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ABSTRACT 


A series of underwater and air explosion investigations was conducted using the 
Arbitrary Lagrangian-Eulerian (ALE) numerical technique. The investigation primarily 
examined the explosive-fluid, fluid-structure, and fluid-air interaction effects, and the 
shock wave pressure propagation through a subjected medium, with the intent of 
verifying and validating the ALE analysis. The research also noted the explosive-air and 
air-structure interaction effects as well as shock wave pressure propagation effects. 
Three-dimensional underwater explosion analyses was conducted using TNT detonations. 
Two-dimensional air explosion analyses were completed using TNT detonations. With 
viable ALE results, underwater and air explosion modeling and simulation could become 
dependable, cost-effective, and time-efficient. 
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I. INTRODUCTION 


A. BACKGROUND 


Air and underwater explosion blast waves have been of great interest to the 
military. An explosion created by a mine produces a pressure pulse or shock wave. 
When a shock wave impacts a structural surface, the shock wave can cause severe, 
negative structural and equipment damage, as well as personnel casualties. As a result, 
military hardware must be shock hardened to ensure combat survivability to both 
personnel and equipment. NAVSEA 0908-LP-000-3010A [Ref. 1] and MEL-S-901D 
[Ref. 2] provide the Navy with guidelines and specifications for shock testing and 
hardening of shipboard equipment. OPNAVINST 9072.2 [Ref. 3] requires a total ship 
system design validation through a series of shock trials. Shock trials are performed on 
the lead ship of a new class of surface ship in order to correct design deficiencies early. 

Unfortunately, shock trials also need extensive planning and coordination. For 
example, planning shock trials for USS Jon Paul Jones (DDG-53) started four years prior 
to the test date. Furthermore, an environmental lawsuit delayed shock trials by three 
months. Shock trials happened in June 1994, but only two of the four planned tests could 
be accomplished due to inclement weather, ship’s post trial delivery date, and deployment 
preparations [Ref. 4]. 

Live Fire Test and Evaluation (LFT&E) assesses a system’s vulnerabilities under 
realistic combat conditions. LFT&E’s objectives are to correct design deficiencies and 
reduce personnel casualties. LFT&E can be expensive due to the destruction and/or 
damage to the ship’s structure and equipment after a shock trial [Ref. 5], Furthermore, 
realistic testing is often sacrificed because of environmental constraints. Thus, realistic 
combat testing of military systems is highly impractical [Ref. 6]. Modeling and 
simulation offers a potential advantage for studying structural failure responses and 
weapons effects under simulated combat test conditions [Ref. 6]. 

Furthermore, land mines represent a serious threat to personnel and motor 
vehicles. Designing a vehicle to withstand the small mine effects is equally important. 
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The value of such a vehicular design can clear lanes through a scatterable minefield. But, 
the need to understand the air blast wave propagation and the resulting mine solid particle 
dynamics is paramount to the vehicular design. Modeling and simulation of the air blast 
wave and the effects of the mine fragmentation can offer such a perspective. 

With recent computer hardware technology advances and increased research in 
numerical simulation of partial differential equations, finite element modeling and 
simulation provides a viable, cost effective alternative to live fire testing under realistic 
combat conditions. Using commercial software packages such as TrueGrid® [Ref. 7] and 
LS-DYNA [Ref. 8] have enabled researchers to produce detailed finite element models in 
a timely manner. "Sufficient fidelity" mesh modeling is required to produce acceptable 
results [Ref. 9]. 

Hydrocodes are defined as tools for the simulation of multi-material, 
compressible, transient continuum mechanics (i.e. mechanical wave propagation through 
multiple fluids and solids) [Ref. 10]. Hydrocodes for mesh descriptions categorically fall 
into two basic arenas: Lagrangian and Eulerian. Each has its distinct advantages and 
disadvantages. Over the past few years, various hybrid hydrocode schemes have been 
developed to handle both Lagrangian and Eulerian. The Arbitrary Lagrangian-Eulerian 
(ALE) numerical formulation was mainly incorporated to study the mitigation effects of 
the explosive-water, explosive-air and structure-fluid/air medium interaction. 


B. SCOPE OF RESEARCH 

This paper investigates the Arbitrary Lagrangian-Eulerian (ALE) formulation 
technique using LS-DYNA in both underwater and air explosion computer modeling and 
simulations. For the underwater explosion (UNDEX) computer simulation modeling, the 
thesis explores the explosive charge modeling effects in shock wave propagation in the 
fluid medium models and the structural response of the structural finite element models. 
Based on the numerical feasibility or confidence of the UNDEX simulation, the thesis 
will then explore the air explosion (AIREX) simulation modeling. The AIREX scenarios 
will investigate the air blast wave propagation and the effects of the blast wave 
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propagation on a varying cylindrical surface thickness. As an overview, Chapter II will 
discuss the important features of the numerical analysis code used in the ALE 
formulation. Chapter III will cover underwater shock principles. Chapters IV will 
describe the modeling setup for the UNDEX scenarios. After gaining numerical 
confidence in the ALE technique, Chapter V will describe the AIREX model setups. 

| 

And, Chapter VI will discuss the results obtained from the conducted simulations. 
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II. NUMERICAL ANALYSIS 

A. ARBITRARY LAGRANGIAN-EULERIAN (ALE) TECHNIQUE 

Explosions involve liquid and gas flow, as well as high-pressure shock waves. A 
Lagrangian finite element mesh in the explosive charge region is not always feasible. 
The surrounding fluid medium elements around the explosive charge deform severely in 
Lagrangian based meshes. Consequently, the time step size per iteration becomes 
extremely small resulting in large computational time [Ref. 11], Furthermore, numerical 
approximation inaccuracies can exist due to mesh distortions [Ref. 12]. 

Eulerian based finite element modeling advance solutions in time on a fixed mesh 
using Navier-Stokes equations. When the solutions are progressed on a fixed mesh, the 
Eulerian hydrocodes avoid mesh distortions as presented in the Lagrangian hydrocodes. 
Additionally, algorithms have been developed to prevent the diffusion between two 
material types at a higher computational expense. Furthermore, solving the Navier- 
Stokes equations (Eulerian) are generally more expensive computationally and 
complicated than the Lagrangian formulation [Ref. 13]. As a result, a hybrid numerical 
formulation technique has been developed which tries to utilize the advantages of both 
the Eulerian and Lagrangian schemes. 

The numerical analysis processor conducted in this study utilizes an ALE finite 
element code. LS-DYNA [Ref. 8] was used for the numerical analysis during this 
investigation. ALE hydrocodes utilize both Lagrangian and Eulerian hydrocodes that 
perform automatic rezoning [Ref. 14]. An ALE hydrocode involves a Lagrangian time 
step followed by a remap or advection phase. The advection phase may pursue one of 
three avenues in which the spatial mesh is (a) not rezoned due to reasonable mesh 
deformation (Lagrangian), (b) rezoned to its original shape due to severe mesh 
deformation (Eulerian), or (c) rezoned to a more suitable form (Lagrangian and Eulerian) 
thus allowing the topology of mesh to remain fixed [Refs. 10 and 14]. It provides 
suitable material models and essential equations of state (EOS) for underwater and air 
explosions. Furthermore, the code provides advection and coupling algorithms in the 
ALE method in order to provide accurate, stable, conservative, and monotonic results. 
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Mass, momentum, and energy transport is systematically computed for all elements in the 
model. Each element’s density, velocity, and energy will be updated. Pressure in each 
element is computed using the updated density and specific internal energies in the 
model’s EOS. 

B. TIME INTEGRATION 

The processor uses the central difference method to advance the mesh position in 
time [Refs. 13 and 14]. The Eulerian time step (At) requires stability and is a function of 
the element’s characteristic length (Ax), material’s speed of sound (c) and particle velocity 
(u) such that 


At< 


Ax 


( 2 . 1 ) 


c + u 

The displacement (x) and velocity (u) vectors are staggered in time to provide a second 
order accurate scheme in time for an "n" iteration step. 


x n + l = x n + u n + ^ ^At n 
U n +1 / 2 _ u n 1/2 + -I a n (At n + At n + 


( 2 . 2 ) 

(2.3) 


F” 

where acceleration vector (a n ) is —. The total nodal force vector and diagonal mass 

M 

matrix are represented by F n and M, respectively. Substituting the acceleration term into 
Equation (2.3) yields: 


u n + l/2 = u n-l/2 4 F 


n 


2M 


(At n + At n + 1 ) 


(2.4) 
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The total nodal force vector consists of an internal nodal force vector (Fj",) and an 


external nodal force vector (F e " t ) such that 


F n = )4 + F e n x t (2.5) 

The internal nodal force is a function is a function of the stress (<r n ) where the stress 
includes the deviatoric sum of the equation of state pressure (-P"I d ) and material strength 
(<j d ) vectors such that 


f^ t = J B ‘o- n dx (2.6) 

V 

cr n = -P n I d + <T d (2.7) 

where B‘ is the strain-displacement matrix, and Id is a principal strain invarient. The 
external nodal force vector (F e " t ) consists of the body forces, boundary forces, non¬ 
reflecting boundary conditions, and contact forces [Refs. 13 and 14]. 

C. MATERIAL MODELS 

This research involves several types of material models. The models incorporated 
the following material models: (a) plastic kinematic/isotropic, (b) null material, and (c) 
high explosive bum. 

The plastic kinematic/isotropic model is used for modeling composites, metals, 
and plastics. This material type effectively measures strain-rate and failure effects and is 
used in beam, shell and solid element modeling [Refs. 8 and 14]. The plastic 
kinematic/isotropic material type is used in modeling the structural metals. 

For the air and water modeling, equations of state (EOS) are employed with the 
null material type. The material strength in the null material element is ignored when the 
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model is combined with a high explosive material type. The null material type is 
effective in modeling fluids and hydrodynamic medium. Additionally, this material type 
can measure failure and thermal effects [Refs. 8 and 14]. The solid element deformation 
is due to a displacement gradient or strain (£). The fluid element deformation is due to a 
velocity gradient or strain rate ( e ). For fluid elements, the deviatoric shear stress (<r d ) is 
proportional to the shear strain rate (£ ) such that 


G d =2H£' ( 2 . 8 ) 

where ji is the fluid vicosity. When a fluid experiences shear stress, the fluid deforms as 
long as a shear stress is applied [Ref. 13]. Hence, air and water are modeled as null 
materials. 

In the high explosive bum model material type, an EOS is used. In this material 
type, bum fractions, F, direct a chemical energy release for detonation simulations. The 
bum fraction is taken as the maximum: 


F = max(^,F^) (2.9) 

where Fi is a function of density, p, the explosive detonation velocity, D, ratio of 

v 

volumetric compression, —, and the Chapman-Jouget pressure, P C j such that 


*L 


pD 


2 ( 


cj 


v 


{ v 0. 


( 2 . 10 ) 


and F 2 is a function of the detonation velocity, D, bum time, tb, current time, t, and 
characteristic length of element, Ax. 
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2 (t - t b )D 
3 Ax 


( 2 . 11 ) 


If the bum fraction, F, exceeds unity, F is reset to one and is held constant [Refs. 8 , 13, 
and 14]. The high explosive pressure, P, in an element is scaled by the bum fraction, F, 
such that: 


P = FPeos( v - E ) < 2 - 12 > 

where P e0 s is the pressure from an EOS based on the relative volume, V, and internal 
energy density per unit initial volume, E [Ref. 14]. Thus, TNT is modeled using a high 
explosive bum material model type. 

D. EQUATIONS OF STATE (EOS) 

An equation relating the pressure, temperature, and specific volume of a 
substance is known as an EOS. Property relations involving other properties of a 
substance at equilibrium states are also known as an equation of state [Ref. 15]. This 
investigation utilized three different EOS in the modeling and simulation. The EOS 
involved were the linear polynomial, Gruneisen, and Jones, Wilkins, and Lee (JWL) 
equations. 

Air is modeled using the linear polynomial EOS. The linear polynomial EOS is 
linear in internal energy per unit initial volume, E. The pressure is given by 

P = Cq + +C^// + C^/^)E (2.13) 


where Co, Ci, C 2 , C 3 , C 4 , C 5 , and C 6 are constants and 



(2.14) 
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where V is the relative volume [Ref. 14]. 

The Gruneisen EOS is used to model the seawater. The Gruneisen EOS 
incorporates a cubic shock velocity-particle velocity, which defines the pressure for a 
compressed material as 


P 0 C 2 JU 




(1 r °] 

a 2 

2 J 

2 


^ + (/ 0 +a//)E (2.15) 


1-(S 1 -1)//-S 2 






f* + l 3 {ju + 1) 2 


where C is the intercept constant of the shock wave velocity (u s )-particle velocity (u p ) 
curve, Yo is the Gruneisen gamma, "a" is the first order volume correction to yo, and Si, 
S 2 , and S 3 are slope coefficients of the u s -u p curve [Ref. 14]. Compression is defined 
using the relative volume, V, as: 


// = “1 (2-16) 

The particle velocity (u p ) is related to the shock wave velocity (u s ) [Refs. 13 and 14] 
through 


U S =C + S l U p +S 2 




v u s; 


U p +S 3 


u, 


'N 3 


V u s J 


u. 


(2.17) 


For TNT, a JWL EOS is used. The JWL EOS defines the pressure as a function 
of the relative volume, V, and initial energy per initial volume, E, such that 
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P = A 


(2.18) 


CO 




R i v ; 


,e R 1 V + B 


( 


co 


K 


RoV 


-RoV 


2 ) 


ClE 

~v 


The parameters (0, A, B, Ri, and R 2 are constants pertaining to the explosive. This EOS 
is well suited because it determines the explosive’s detonation pressure in applications 
involving structural metal accelerations [Ref. 17]. 

E. ADVECTION ALGORITHM 

A good advection (remap) step involves accuracy, stability, conservation, and 
monotonicity [Ref. 14]. A monotonic, second order accurate "Van-Leer and Half-Shift 
Index" advection scheme is used for material transport. In mass advection, a new 
Eulerian density for each element is solved using the following scheme [Ref. 13]: 


Pe V e = PL V L+ 2 W/j (2 ' 19) 

faces 


where 

p e - Eulerian element density 

V e - Eulerian element volume 

p L - Lagrangian element density 

V L - Lagrangian element volume 


- Density of adjacent Lagrangian element j 
4 >j - Volume flux through adjacent element j 


Internal energy advection is expressed in terms of a variable, s, where 


E 

s = — 
V 


( 2 . 20 ) 
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( 



Or, mass advection can be expressed in terms of a variable, s, such that, 


s = P 


( 2 . 21 ) 


Thus, equation (2.19) is now written as a general advection equation of s e V e where 


s e V e“ s L V L + ^ s i^j 
faces 


( 2 . 22 ) 


In a one-dimensional coordinate system, the change in momentum in the x- 
direction can be solved in terms of the u-velocity. The Lagrangian element u-velocity 

(uf! em ) can be expressed as 


?> em =l | u i 


U L 


«j=r L 


(2.23) 


plpm 

Using the velocity result from Equation 2.23, the variable x-momentum (s^ ) is 


p elem .elem „elem 
S L =PL u L 


The updated or new element centered x-momentum flux becomes 


(2.24) 


s e V e = s L V L + AM ) 


(2.25) 


where AM X is defined as the change in x-momentum resulting from advection: 
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(2.26) 



A similar analysis is applied to Equations (2.19) through (2.26) when dealing with a 
three-dimensional model. 

Using the analysis in Equations (2.23) - (2.26) yields a first order accurate 
scheme. The Van-Leer scheme yields a second order accurate advection scheme where 
cl cm 

the variable, Sl , is modified (Figure 1). Instead of using Equation (2.24), the Van- 

cl cm 

Leer scheme utilizes Equation (2.27) to calculate : 

SL em = / , O + |( x l- x 0)^( x l) < 2 - 27 > 

F. COUPLING ALGORITHM 

In Figure 2, the structure-medium coupling interface designates a structural 
surface, typically shell elements, as a "slave" material and the fluid as a "master" material 
[Refs. 13 and 14]. This model uses a “penalty coupling” factor in the Lagrangian- 
Eulerian coupling. No special pre-processing manipulation of the intersecting surfaces 
between the master and slave materials is required. 

The penalty coupling factor tracks the relative displacement (d) between a 
Lagrangian node (structure i.e. slave material) and the Eulerian fluid (master) material 
location (Figure 3). Each slave node is checked for penetration through the master 
surface. If no slave node penetration occurs, nothing is done. If slave node penetration 
into the master surface occurs, an interface force (F) is distributed to the Eulerian fluid 
nodes. The interface force's magnitude is proportional to the amount of penetration 
occurring such that: 


F = k i -d 


(2.28) 
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where kj is a stiffness factor based on the master and slave nodes mass model properties. 

As the interface force (F) in Equation (2.28) is solved for each time integration, F 
is considered as one of the external body forces in Equation (2.5). Thus, a total nodal 
force (F n ) can be solved in each time integration resulting in structural accelerations, 
velocities, and displacements. 
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Figure 1. Second Order Van Leer Algorithm (From Ref. 13) 
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Figure 2. Lagrangian/Eulerian Coupling (After Ref. 13) 
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timet time t+At 



Figure 3. Penalty Coupling (From Ref. 13) 
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III. THEORY 


A. UNDERWATER SHOCK PHENOMENA 


An underwater explosion is a complex sequence of events. The explosion starts 
with a high explosive detonation. The explosion consists of a chemical reaction in an 
initial unstable substance. The initial chemical reaction converts the unstable material 
into a stable gas at a very high temperature and pressure (on the order of 3000° C and 
50000 atm.) [Ref. 18]. As the initial pressure wave within the explosive disseminates, 
additional chemical reactions follow thus creating additional pressure waves. The 
pressure wave velocity steadily increases within the explosive until velocity exceeds the 
speed of sound in the explosive. Consequently, the explosive produces a shock wave. 
The combination of high temperature and pressure describes the detonation as a self- 
exerting process where the temperature and pressure profiles exist behind the shock 
wave. The high temperature and pressure gas expands in an outward motion into the 
surrounding fluid medium [Ref. 19]. 

The high-pressure disturbance traveling outward radially at the speed of sound is 
referred to as the shock wave. As the shock wave expands outward, the explosion also 
induces an outward velocity to the surrounding fluid medium as described in Figure 4. 
The initial high-pressure is greater than the opposing atmospheric and hydrostatic 
pressure. Thus, the initial pressure is compressive in nature. At detonation, the pressure 
rise produces a steep fronted discontinuous wave decaying exponentially with respect to 
time as shown in Figure 5. The pressure disturbance lasts only a few milliseconds. 
Although the shock wave character near the explosive charge is unknown, the shock 
wave propagates roughly at several (3—5) times the speed of sound in water, 
approximately 5,000 ft/sec, which then falls rapidly to acoustic velocity as it travels 
outward in the water [Ref. 18]. 






High Pressure Gas Bubble 
| Charge Diameter 

ShockWave Front 


Figure 4. Gas Bubble and Shock Wave from an Underwater Explosion 

[From Ref. 19] 


34,000 psi 



3,400 psi 




495’ 500' 


Figure 5. Shock Wave Profiles From a 300 lb. TNT Charge [ From Ref. 18] 




Furthermore, the shock wave pressure is proportional to the inverse of the 
distance from the charge to the target, 1/d. Additionally, the shock wave profile 
gradually broadens as it spreads out [Ref. 18]. Empirical equations have been formulated 
to describe the shock wave profile. These relations calculate of the shock wave pressure 
profile, P(t), the maximum shock front pressure, P mM , the shock wave decay constant 0, 
the gas bubble period (T), and the maximum gas bubble radius (A max ). 


-Mi) 

p ( t ) = P max e e (psi); t> tj (3.1) 


*max *4 


( n A i 

W3 
R 

V J 


(psi) 


(3.2) 


1 

e =K 2 w3 


( n A 2 

W3 
R 

V ) 


(msec) 


(3.3) 


T = K 


1 

W3 
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(D + 33) 6 


(sec) 


(3.4) 
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^max ^6 


1 

W3 


(ft) 


(3.5) 


1 

(D + 33) 3 


where: 

W = Charge weight (lbf) 

R = Standoff distance (ft) 

D = Charge depth (ft) 

tj = arrival time of shock wave (msec) 

t = time of interest (msec) 

Ki, K 2 , K 5 , K 6 , Ai, A 2 = Explosive Shock wave parameters 


Equation (3.1) is good only for pressures greater than one third of P max (Equation 3.2). 
Through calculation, it can be determined that Pmax decreases by approximately one-third 
after one decay constant ( 0 ). 

In the underwater explosion process, the gas bubble oscillation produces 
subsequent pressure waves or bubble pulses. The peak pressure of the first bubble pulse 
is approximately 10 - 20 % of the shock wave, but is of greater duration making the area 
under both pressure curves similar [Ref. 19]. The bubble expands until dynamic 
equilibrium is reached . The bubble then contracts until dynamic equilibrium is again 
reached, followed by another expansion. This oscillation sequence continues until the 
energy of the reaction is dissipated or the bubble reaches the free surface or impacts the 
target. 

Depending on the charge location relative to the surface and the bottom, other 
issues have an effect on underwater shock. Bottom reflection waves involve bouncing 
shock wave off the sea bottom. Bottom reflection waves are considered compressive 
waves in nature. Refraction waves encompass shock waves traveling through the sea 
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bottom before emerging again in the water. Refraction waves are also considered 
compressive waves. In deep water, these two phenomena are not usually an issue for 
surface vessels. 

Free surface reflection is a very important effect, however. Reflection or 
rarefaction waves are tensile, as opposed to the other compressive wave effects. 
Rarefaction waves are produced from the shock wave reflecting from the free surface. 
Consequently, rarefaction waves contribute to the bulk cavitation phenomenon. 


B. CAVITATION 


During an UNDEX, two types of cavitation can occur, namely local and bulk 
cavitation. Local cavitation occurs at the fluid-structure interface, while bulk cavitation 
occurs near the free surface and covers a relatively large area. 


1. Local Cavitation 

Taylor flat plate theory, the simplest case of fluid-structure interaction, depicts 
how local cavitation occurs. In this case, an infinite, air-backed plate is acted upon by an 
incident plane shock wave as shown in Figure 6 . 

Once the shock wave strikes the plate, a reflected shock wave leaves the plate. 
Using Newton’s second law of motion for the equation of motion yield: 

m^ = Pi+P 2 (3.6) 

dt 1 z 
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Figure 6. Taylor Plate Subjected to a Plane Wave [From Ref. 19] 

where m is the plate mass per unit area, u is the plate velocity after being subjected to the 
shock wave, Pj(t) is the incident wave pressure and P 2 (t) is the reflected, or scattered, 
wave pressure. The fluid particle’s velocities behind the incident and reflected shock 
waves are defined as uj(t) and u 2 (t). The plate velocity, u(t), is defined by Equation (3.7), 


u(t) = U^t) - u 2 (t) 


(3.7) 


For the one dimensional plane wave, the wave equation is P = pcu. The incident 
and reflected shock wave pressures are shown in Equations (3.8) and (3.9) respectively: 


P^pcuj 


(3.8) 


F 2 = pcu 2 


(3.9) 


where p is the fluid density and c is water’s acoustic velocity. Substituting Equations 
(3.8) and (3.9) into the plate velocity Equation (3.7) yields: 
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u = 


P l-P2 


(3.10) 


pc 


Using Equation (3.1), the incident shock wave pressure, Pi(t), now becomes: 


-t 

P 1 = P max e ^ (3-11) 

where time, t, is measured after the shock wave arrives at the target. Using Equation 
(3.10), the reflected shock pressure ,P 2 (t), gives: 


P 2 =P 1 -pcu (3.12) 


-t 

P 2 = P ma x * e -PC" (313) 


Combining the equation of motion, Equation (3.6), with Equations (3.11) and (3.13): 


du 

m-1-pcu 

dt 



(3.14) 


du OT> 

m d! +pcu=2Pmax ' 


-t 

7e 


(3.15) 
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where Equation (3.15) is a first order, linear differential equation. Solving the 
differential equation in Equation (3.15) produces a solution, u(t): 




r-A) 


(—) 

^max^ 

e 

l e J 

— e 

< 0 j 


m(l-0) 


with p = pC0/m and t>0. The total pressure that impinges on the plate is defined as: 


Pi + P 2 P max 


1-/5 


-t 

7e 


~(5t 

1 -fi 


(3.17) 


As the value of p becomes larger, as in the case of a lightweight plate, the total pressure 
will become negative at a very early time. However, since water cannot support tension, 
negative pressure cannot exist. Therefore, as the water pressure reduces to vapor 
pressure at the surface of the plate, local cavitation occurs. At this point, the pressure in 
front of the plate has been cut off and the plate has reached its maximum velocity [Ref. 
19]. 

A ship’s hull can be easily generalized as a Taylor flat plate. Local cavitation is 
likely to occur along the hull where the pressure pulse from the UNDEX impinges with 
sufficient force and the hull plating p value is large enough to make the net pressure 
negative. 
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2 . 


Bulk Cavitation 


As discussed previously, the incident shock wave is compressive in nature. A 
rarefaction wave (tensile wave) is created when the shock wave is reflected from the free 
surface. Since water cannot sustain a significant amount of tension, cavitation occurs 
when the fluid pressure drops to zero or below. Upon cavitation, the water pressure rises 
to the vapor pressure of water, approximately 0.3 psi. This cavitated area produced by 
the rarefaction wave is known as the bulk cavitation zone. The bulk cavitation region 
contains an upper and lower boundary. The bulk cavitation zone relies on the charge 
size, type, and depth [Refs. 20 - 22]. 

Figure 7 shows a typical bulk cavitation zone. The cavitation zone is symmetric 
about the y-axis in the figure; typically only one-half is shown due to the symmetry about 
the y-axis. At the time of cavitation, the water particles’ velocities behind the shock wave 
front depend on their location relative to the charge and the free surface. For example, 
the water particles’ velocities near the free surface will have a primarily vertical velocity 
at cavitation. As the reflected wave passes, gravity and atmospheric pressure will 
primarily act on the water particles’ velocities. 

The upper cavitation boundary consists of a set of points where the rarefaction 
wave passes and reduces the absolute pressure to zero or a negative value. This region 
will remain cavitated as long as the pressure remains below the vapor pressure. The total 
or absolute pressure, which determines the upper boundary, is composed of atmospheric 
pressure, hydrostatic pressure, incident shock wave pressure, and rarefaction wave 
pressure. 
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Figure 7. Bulk Cavitation Zone [From Ref. 19] 
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The lower cavitation boundary is determined by equating the decay rate of the 
breaking pressure to the decay rate of the total absolute pressure. The breaking pressure 
is the rarefaction wave pressure that reduces a particular location of a fluid to the point of 
cavitation pressure, or zero psi. 

The upper and lower cavitation boundaries are calculated using Equations (3.18) 
and (3.19), respectively [Ref. 21]. Any point which satisfies F(x,y) and G(x,y) = 0 
determines the bulk cavitation boundary. 
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The variables in Equations (3.18) and (3.19) are: 


x, y = horizontal range and vertical depth of the point 

ri = standoff distance from the charge to the point 

r 2 = standoff distance from the image charge to the point 

c = acoustic velocity in the water 

D = charge depth 

0 = decay constant 

y = weight density of water 

P A = atmospheric pressure 

W = charge weight 

0 = decay constant 

Ki, A] = shock wave parameters 


Figure 8 shows the charge geometry for Equation (3.18) and (3.19). 

Appendix A provides a MATLAB m-file [Ref. 23] calculating and plotting the 
bulk cavitation zone for a user supplied charge weight (HBX-1 and TNT) and depth by 
solving Equations (3.18) and (3.19). Figures 9 and 10 provide an example of bulk 
cavitation curves generated using the program for two different charge depths. 
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Cavitation zone for 60 lb TNT Charge 



Figure 9. Bulk Cavitation Zone for 60 lbf TNT Charge at Depth of 25 ft 
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IV. UNDERWATER EXPLOSION 


Modeling and simulation involved a pre-processor, processor, and post-processor. 
A flow chart of the model building and testing procedure is shown below in Figure 11. In 
this investigation, the models were generated or preprocessed using TrueGrid. Next, 
LS-DYNA was used for the analysis and solution (processor). Finally, results were 
displayed utilizing LS-POST and MATLAB programs. 



Post-Processor 

LS-POST 

MATLAB 


Figure 11. Flow Chart Model Construction and Simulation 
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A. 


MODEL CONSTRUCTION AND DIMENSIONS 


TrueGrid was used to model several underwater explosion scenarios. The TNT 
explosive was modeled first, since the TNT’s geometry dictated the entire model mesh. 
Next, the surrounding fluid and air meshes were generated. And finally, a floating 
structural model, a barge, was modeled as an "independent system" within the fluid mesh. 
Appendix B provides a general listing of TrueGrid commands [Ref. 7]. 

1. Explosive Model 

The underwater explosion scenarios started with explosive modeling. Two 
different models were used in this study. The two models involved the use of a TNT 
explosive. The models were differed by the explosive geometry, namely a block charge 
and a “stair-cased” or “stair” charge. Since high explosives react rapidly to provide 
gaseous products at high temperatures and pressures from an initial volume, the explosive 
was modeled from an initial charge density (p) and charge weight (W). A charge weight 
of 60 lbf or 266.9 N was used in the TNT models. The explosive’s initial densities were 
obtained [Refs. 17 and 24], Knowing 


W = mg (4.1) 

m 

V = — (4.2) 

P 

where g is the gravitational constant, m is the explosive mass, V is the explosive’s initial 
volume. So, the TNT's initial volume was computed. The TNT explosive was initially 
modeled as a rectangular box (Figure 12). Appendix C provides the TrueGrid commands 
for the TNT rectangular explosive block case. 

The entire mesh for the rectangular explosive block case would prove to be too 
coarse during the simulation runs, so the entire mesh was refined through remodeling the 
TNT explosive. A stair-cased shape TNT explosive model was generated to simulate a 
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quarter of a sphere. In order to accomplish this, the original volume, V, was broken up 
into fourteen separate sub-volumes, V*, such that 


* 

V 


_v 

14 


(4.3) 


These sub-volumes were then arranged in a “stair” manner to simulate a one- 
fourth of a spherical mine (Figures 13 through 15). Furthermore, the “stair” element 
length size, l s , was calculated from V* such that: 



Since stair element length size is less than the rectangular block element length size, the 
overall mesh modeling is much more refined with the stair element length size. 
Appendix D provides the TrueGrid commands for stair TNT explosive. Table 1 provides 
the charge modeling characteristics for both the rectangular block and stair geometries. 

The TNT rectangular block explosive charge consisted of 16 nodes, which 
composed 3 eight-noded solid elements. For stair explosive charge, the explosive model 
consisted of 50 nodes, which composed 14 eight-noded solid elements. LSDYNA's 
MAT_HIGH_EXPLOSIVE_BURN material type was used to model the TNT. Table 2 
illustrates TNT's material properties [Ref. 17]. TNT material type also requires an EOS. 
Equation (2.18) is the JWL EOS for high explosives. Table 3 provides TNT's parameters 
for the JWL EOS [Ref. 17]. 

Appendices E and F are the LSDYNA input decks for the rectangular box and 
stair charges respectively. The explosive required axis-symmetric boundaries along two 
faces of the explosive. By placing an axis-symmetric boundary along the charge's X-Y 
plane, the nodes along the boundary were translationally and rotationally constrained in 
the Z direction. Additionally, by placing an axis-symmetric boundary along the X-Z 
plane, the nodes along the boundary were translationally and rotationally constrained in 
the Y direction. Figures 12 and 14 illustrate the axis-symmetric boundaries. 
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Charge 

Type 

Charge 

Weight 

(N) 

Charge 
Density, p 
(g/cm 3 ) 
[Ref. 17] 

Initial 

Volume, 

V 

(cm 3 ) 

Charge 

Geometry 

Block 
Element 
length, l b 
(cm) 

Stair 
Element 
length, l s 
(cm) 

TNT 

266.9 

1.63 

16690 

Block 

17.45/ 

18.27 

N/A 

TNT 

266.9 

1.63 

16690 

Stair 

N/A 

10.60 


Table 1. TNT Charge Modeling Characteristics 


Charge 

Detonation Velocity, D 
(cm/ps) 

Chapman-J ouget 

Pressure, P C j (Mbar) 

TNT 

0.693 

0.210 


Table 2. TNT Charge Characteristics (Ref. 17) 


Charge 

A 

(Mbar) 

B 

(Mbar) 

Ri 

r 2 

CO 

Eo 

Mbar- 

cm 3 /cm 3 

Vo 

(cm 3 / 

cm 3 ) 

TNT 

3.71 

3.23 x 
10' 2 

4.15 

0.95 

0.30 

4.30 x 
10‘ 2 

1.0 


Table 3. TNT Equation of State Parameters (Ref. 17) 
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2 . 


Fluid Model 


Upon completion of the explosive model, the fluid mesh was designed next. 
Appendices C and D provide the TrueGrid commands for the fluid mesh. The fluid mesh 
was primarily designed around the respective explosive meshes. Upon initial fluid mesh 
construction, several fluid solid element blocks were deleted to fit the explosive mesh 
adjacently to the fluid. Furthermore, additional fluid solid element blocks were deleted 
on the water and air surface boundary. These fluid blocks were removed to make room 
for the barge placement and air solid element blocks within the barge. Fluid nodes would 
serve as an “interface” which prevents the fluid nodes from merging with the structural 
surface (Figure 18). As with the explosive charge, axis-symmetric boundaries were 
required along the fluid's X-Y and X-Z planes. The nodes along the X-Y and X-Z planes 
were translationally and rotationally constrained in the Z and Y directions, respectively. 
Along three of the fluid model sides, non-reflective boundaries were placed to prevent 
unwanted reflected "noise" from those sides during the simulation runs. Figures 16 and 
17 illustrate the fluid mesh’s boundary constraints. Along the sixth fluid side, namely the 
fluid-air medium interaction, no constraints were placed. 

LS-DYNA's NULL_MAT is used to model the fluid material. A water density of 
1.0 g/cm 3 was used. The Gruneisen equation of state (Equation 2.15) is also incorporated 
to model the water. The shock Hugoniot parameters for water were obtained from the 
HULL code [Ref. 16] where 

C = 0.148 cm/|xs 

S, = 1.75 

s 2 s 3 0 

Yo = 0.4934 

The overall dimensions for the fluid mesh are 1036.32 cm x 731.52 cm x 670.56 
cm. The fluid mesh density was primarily based on the explosive element length size for 
each scenario. In the TNT rectangular block explosive model, the fluid model possessed 
100474 nodes, which composed 93961 eight-noded solid elements. In the stair model, 
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the fluid model consisted of 1573224 nodes, which composed 148518 eight-noded 
elements. 


3. Air Model 

After the fluid mesh was completed, the air mesh was constructed. Appendices C 
and D show the TrueGrid commands for the air meshes. After initial construction of the 
air mesh, air solid element blocks were deleted to fit the fluid mesh adjacently with the 
air mesh. Additionally, an axis symmetric boundary was placed along the X-Y plane 
requiring the nodes along the X-Y plane to be translationally and rotationally constrained 
along the Z direction. Along three of the air model sides, non-reflective boundaries were 
placed. Figures 20 and 21 describe the air mesh’s boundary conditions. On top of the air 
model and on the fluid-air medium interaction surface, no boundary constraints were 
placed. 

LS-DYNA's NULL_MAT is used to model the fluid material. An air density of 
1.28*1 O' 3 g/cm 3 was used. The linear polynomial EOS (Equation 2.13) is used to model 
the air. The air was modeled as a near perfect gas such that the constants for Equation 
(2.13) were 


Co = 0 

Cl = 0 

C 2 = 0 

C 3 = 0 

C 4 = 0.4 

C 5 = 0.4 

C 6 = 0 

The overall general dimensions for the air mesh are 1036.32 cm x 182.88 cm x 
670.56 cm. The air mesh density was also based on the explosive element length size for 
each scenario. In the rectangular block explosive case, the air model possessed 26000 
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nodes, which composed 22620 eight-noded solid elements. For stair explosive case, the 
air model consisted of 41250 nodes, which composed 36664 eight-noded elements. 

4. Structural Model 

Finally, a ship-like box model, barge, was constructed for ship shock simulation. 
The barge’s overall dimensions are 426.72 cm x 60.96 cm x 60.96 cm (Figure 24). The 
barge possessed two athwartships bulkheads (Figure 25). The barge was weighted with 
three lumped masses (1.163 x 10 s g) evenly spaced along the barge’s keel ensuring the 
center of gravity remained along the centerline (Figure 26). The barge’s shell plating 
consisted of 0.635 cm thick steel, 0.2% C hardened, having a weight density 7.87 g/cm , 
a Young’s modulus of 210 GPa, a Poisson’s ratio of 0.3, and a yield strength of 430 MPa 
[Ref. 24]. Furthermore, an “interface” was developed along the barge’s outer structural 
nodes in order to prevent the nodes from merging with the fluid nodes. The overall 
structural finite element mesh consisted of 443 nodes, which composed 432 four-noded 
thin shell elements. Appendices C and D provide the TrueGrid commands for the 
structural model inputs. 


B. ANALYSIS AND SOLUTION 

1. Analysis Code Description 

TrueGrid translates the generated finite element model into an LS-DYNA 
keyword format in order to perform the numerical analysis. LS-DYNA's arbitrary 
Lagrangian-Eulerian (ALE) numerical technique was incorporated in the shock analysis. 
Utilizing LS-DYNA's *CONSTRAINED_LAGRANGE_IN_SOLID command couples 
the independent Lagrangian mesh (structure shells) with the independent Eulerian mesh 
(fluid and air solid elements). A penalty-coupling factor of 0.1 is used for the simulation 
runs. Chapter II discusses the analysis and solution for the finite element mesh. The two 
different models were set to run for 10 msec. The finite element models used a time step 
scale factor of 0.67 because high explosives were used [Ref. 25]. Appendices E and F 
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provide the LS-DYNA input decks for the rectangular block and stair explosive cases, 
respectively. Appendix G provides useful LS-DYNA commands. Lastly, the units used 
in LS-DYNA were in centimeters, grams, and microseconds. 


2. Test Description 

Two different charge models, rectangular block and stair, were used in the shock 
simulation runs for this study. The explosive consisted of a 266.9 N charge. The attack 
geometry placed the charge offset from the side of the barge model by 304.8 cm. The 
charge depth was 731.52 cm, with a standoff distance of 792.48 cm. Figure 27 shows the 
attack geometry for this offset charge scenario. Using Appendix A, the bulk cavitation 
zone was computed and is shown in Figure 28. 

C. POST-PROCESSING 

The solution data output is placed into two formats for analysis: binary and 
ASCn. The binary data files created by the LS-DYNA simulations contain the model’s 
finite element response information. The binary data files recorded data every 200 psec. 
Thus, the simulation runs produced 50 states during the 10000 psec time interval. LS- 
POST [Ref. 26] is used for three-dimensional response visualization. The powerful 
animation and image generation features produced by LS-POST include displacement, 
velocity, acceleration, and element pressure data display. Furthermore, LS-POST 
enabled the user to observe the shock wave propagation through the fluid or air medium. 
LS-POST also has the capability of extracting ASCII solution data and writing it to a 
separate ASCII file for later evaluation. The barge’s nodal velocities were of primary 
importance. Eight specific nodes on the barge were recorded using the 
DATABASE_HISTORY_NODE command in LS-DYNA. This feature was crucial in 
obtaining the structural velocity and fluid pressure response and converting the data into 
a usable form for data interpretation. Appendix H provides some useful LS-POST 
commands for model post-processing. The extracted ASCII data by the LS-POST post- 
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processor was plotted and manipulated using MATLAB [Ref. 23]. The graphical output 
was converted from metric to English units while using MATLAB. 
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Figure 14. Stair-Cased Explosive Model for TNT (Side View) 







Figure 15. Stair-Cased Explosive Model for TNT (Angle View) 
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Figure 17. Fluid Mesh (Side View) 
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Figure 18. Three-Dimensional Fluid Mesh 
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Figure 20. Air Model (Top View) 
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Figure 21. Air Model (Side View) 
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Figure 22. Three-Dimensional Air Mesh 
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Figure 23. Three-Dimensional Air Mesh Parameters 
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Figure 24. Three-Dimensional Barge 








Figure 25. Three-Dimensional Barge (Bulkheads) 
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Figure 27. Offset Charge Test Geometry 
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Cavitation zone for 60 lb TNT Charge 



Figure 28. Bulk Cavitation Zones for 60-lbf TNT Charge at 24 ft 
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V. AIR EXPLOSION 


A. MODEL CONSTRUCTION AND DIMENSIONS 

The same process and technique applied to the underwater explosion simulations 
are applied to the air explosion simulations. In this two-dimensional air explosion 
simulation, one-fourth of a cylindrical TNT explosive is modeled. After modeling the 
explosive, the surrounding cylindrical air mesh is constructed, and finally a structural 
shell model is fabricated. 

1. Explosive Model 

As previously mentioned in Chapter IV, a similar methodology is approached to 
determine the explosive’s modeling dimensions. In this two-dimensional air explosion, 
the simulation runs involved a 0.5 lbf or 2.22 N TNT cylindrical charge. Utilizing 
Equations (4.1) and (4.2), the TNT’s initial volume, V, was calculated. One-fourth of the 
cylindrical explosive was modeled in order to save computational time and cost. The 
explosive's radial element length, r, was determined from the initial volume, V. 
Additionally, the radial element length was kept the same as the cylindrical height. So 
Equation (5.1) provides the radial element length such that 



Therefore, in this case the explosive’s radial and height dimensions were 3.54 cm 
(Figures 29 and 30). See Appendix I for TrueGrid command inputs. 

The TNT explosive model consisted of 124 nodes, which composed 61 solid 
elements. LS-DYNA's MAT_HIGH_EXPLOSIVE_BURN and JWL EOS were used in 
conjunction for modeling the explosive material. TNT’s material properties are listed in 
Tables 1 through 2. Figure 30 shows two boundary constraints at the extremes of the 
explosive's radial arm faces in order to maintain the physics of the scenario. 
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2. Air Model 


After the fabricated explosive model, the air model was built next. The air 
model's dimensions were 311.52 cm from the explosive’s center and a depth of 3.54 cm 
(Figure 31). LS-DYNA’s NULL-MAT and the linear polynomial equation of state 
(Equation 2.13) were used to model the air. An air density of 1.28*10' 3 g/cm 3 was used. 
Since perfect gas law assumptions were used, the linear polynomial equations of state 
constants were Co = Ci = C 2 = C 3 = C 6 = 0 and C4 = C 5 = 0.4. As in the explosive model, 
the two boundary constraints were placed at the extremes of the radial arms. A non¬ 
reflecting boundary was placed on the outer-most radial curvature. Figure 32 illustrates 
the boundary conditions of the air model. The air model possessed 10952 nodes, which 
composed 5307 eight-noded solid elements. Appendix I furnishes the TrueGrid 
command inputs. 

3. Structural Model 

A cylindrical structural plate was constructed for shock simulation testing. 
Three different material sizes or thicknesses were modeled. The same material type used 
in Chapter IV is studied in this air explosion. Table 4 provides the steel material's 
characteristics [Ref. 24], The three structural thicknesses investigated were 0.635 cm, 
1.27 cm and 2.54 cm (Figure 33). An air gap of 2.54 cm existed between the explosive 
and cylindrical steel plate (Figure 34). The plate was not weighted and none of the steel's 
nodes were constrained. The overall structural finite element mesh consisted of 60 
nodes, which composed 29 four-noded thin shell elements. Appendix I provides the 
TrueGrid command inputs. 
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Material 

Type 

Density 

(g/cm 3 ) 

Young’s 

Modulus 

(GPa) 

Poisson’s 

Ratio 

Yield 

Strength 

(MPa) 

Ultimate 

Strength 

(MPa) 

Structural 

Steel 

7.87 

250 

0.30 

430 

620 


Table 4. Material Property Characteristics 

B. ANALYSIS AND SOLUTION 

1. Analysis Code Description 

The finite element model was translated into LS-DYNA keyword format in order 
to perform the numerical analysis. LS-DYNA’s arbitrary Lagrangian-Eulerian (ALE) 
numerical technique was instituted in the air shock analysis. Utilizing LS-DYNA’s 
*CONSTRAINED_LAGRANGE_IN_SOLID command coupled the Lagrangian mesh 
(structure shells) with the Eulerian mesh (air solid elements). The penalty-coupling 
factor remained at 0.10 for all simulation runs. Chapter II discusses the analysis and 
solution for the finite element mesh. The three simulation runs were set for 4 msec. The 
finite element models used a time step scale factor of 0.67 or less because high explosives 
were used [Ref. 25]. Appendices J through L provide the LS-DYNA input decks for each 
of the simulation runs. Appendix K shows useful LS-DYNA commands. The units of 
the data obtained from LS-DYNA are in centimeters, grams, and microseconds. 


2. Test Description 

Three different structural thicknesses were used in the air shock simulation runs. 
The explosive consisted of a 2.22 N (0.5 Ibf) charge. The attack geometry places an air 
gap of 2.54 cm exists between the charge and the cylindrical structural plate. Figure 34 
provides the attack geometry for the simulations with varying shell thicknesses. 
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c. 


POST-PROCESSING 


The solution data output was placed into two formats for analysis: binary and 
ASCH. The binary data files created by the LS-DYNA simulations contained the 
model’s finite element response information. The binary data files recorded every 50 
(isec. LS-POST was used for three-dimensional response visualization in displacements, 
velocities, accelerations and pressures [Ref. 26]. Furthermore, LS-POST enabled the 
user to observe and record the shock wave propagation the air medium. Animation 
sequences of the structural response and shock wave pressure propagation were created 
using the movie dialog and saving as an AVI (rle) format. LS-POST also extracted 
ASCH solution data and wrote the data to a separate ASCII file for later evaluation. 
Appendix L provides some useful LS-POST commands for model post-processing. The 
extracted ASCII data by the LS-POST post-processor was plotted and manipulated using 
MATLAB [Ref. 23]. The graphical output was converted from metric to English units 
while using MATLAB. 
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Figure 29. Explosive Model (Angle View) 
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Figure 31. Air Model (Top View) 
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Figure 34. TNT Explosive and Steel Shell Plating 
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VI. SHOCK SIMULATION RESULTS 


All simulations were run on an SGI Octane with a 195 MHz processor, 1.344 
Gigabytes of RAM, and 23 Gigabytes of hard drive storage capacity. LS-DYNA version 
960 (alpha) revision 1.174 was the numerical processing simulation codes used for the 
three-dimensional underwater and two-dimensional air explosions. 


A. UNDERWATER EXPLOSION 


Two simulations were observed for comparison, namely, the difference in 
explosive charge modeling. One simulation involved the use of a rectangular box 
explosive model. The other simulation encompassed the use of a stair-cased shape 
explosive model. In the fluid medium, sixteen solid fluid elements were used for 
comparison in the simulation, specifically five fluid elements on three different faces of 
the fluid mesh (X-Y axis symmetric boundary, non-reflective boundary 1, and non- 
reflective boundary 2) and one fluid element adjacent to the barge. Figure 16 illustrates 
the locations of the three faces or planes in the fluid mesh. The pressure profiles were 
analyzed and compared to the empirical formulations discussed in Chapter EL All of the 
elements exhibited a typical exponential decay waveform. When the shock pressure 
wave arrived at an element at a certain time, the pressure ascended quickly to a maximum 
pressure, P max . Then, the pressure waveform exponentially decayed until the end of the 
simulation run. Furthermore, eight nodal points on the barge (four nodes) and the fluid 
mesh (four nodes) were used for comparison. Due to the bulk cavitation, the vertical 
velocity response, kickoff velocity, was analyzed at the structural and corresponding fluid 
nodes. Figures 35 and 36 illustrate the barge nodal locations for both the rectangular box 
and stair shape cases. The corresponding fluid nodal locations are shown in Figures 37 





1 . 


Rectangular Block Explosive Simulation 


Exploring the water’s pressure response within the fluid model, sixteen solid fluid 

elements were used for comparison. Along the X-Y plane, nodes 1128, 1139, 46247, 
46256, and 46256 were used for comparison to the underwater explosion empirical 
equations. Along the non-reflective boundary (NRB) 1 plane, nodes 67155, 69178, 
79616, 79786, and 89068 were compared to the empirical formulations. And along the 
NRB 2 plane, nodes 25576, 43680,46087,47694, and 47706 were used for comparison. 
The last fluid element involved 46968 was located underneath the barge’s centerline 
approximately. Figures 39 through 54 present the numerical and empirical pressure 
results. The empirical and numerical pressure profiles exhibit similar contour agreement. 
Tables 5 through 7 compare the results in maximum peak pressure, P max , between the 
empirical and numerical formulations. In comparison to the empirical P max , the 
numerical results show relatively good agreement at close locations and less agreement at 
farther locations. Three things may have a numerical effect on the water pressure 
profiles: (1) the explosive charge may need more refinement (as will be seen in the next 
section), (2) the time scale factor may require further reduction in order to capture the 
sharp rise in pressure in the short amount of time, or (3) the shock Hugoniot parameters 
for water may require further investigation. An analysis may be required to compare the 
linear HULL code constants in the Gruneisen EOS with other shock Hugoniot parameter 
constants such as Steinberg’s non-linear parameters [Ref. 16]. 

Cavitation can also be observed in several nodes, namely nodes 46264, 47706, 
and 46968, which were at depths of 2.15 ft, 2.73 ft, and 1.0 ft respectively. The 
mentioned nodes at their respective depths were well within the cavitation zone and were 
in agreement with Figure 28. The other nodes were outside of the cavitational region. 
As expected, the water particle velocities at the fluid-air boundary behaved principally in 
a vertical manner due to the bulk cavitation. Additionally, the air pressure near the 
vicinity of the fluid-air boundary experienced little or no relative change in pressure as 
predicted. 


72 



Since the nodes in the fluid mesh exist independently of the structural mesh, a 
correspondence between the structural nodes and fluid nodes must be established. The 
structural nodes 124024, 124134, 124173, and 124190 correspond to the following fluid 
nodes, respectively, 16862, 39870, 39868, and 40492. In some instances, the structural 
node location did not match up exactly with the fluid node location, so the closest fluid 
node was chosen. Figures 55 through 58 illustrate the rapidly increasing structural and 
water particle vertical velocities (approximately 5 to 6 msec) to a peak value and then a 
rapid decrease and increase in a series of steps until the response settles down. The 
response does not settle out at a value of zero due to rigid body motion of the structure. 
The rigid body motion is due to the fact that the barge exists independently from the fluid 
mesh. The incident shock wave impacts the structure with a very high pressure 
(approximately 2500-psi at 5 msec) from charge detonation and forces the structure 
rapidly upward. The structure is then quickly pulled downward as the shock wave 
reaches the free surface and a rarefaction wave (tensile) is generated. This wave causes 
the pressure to decrease rapidly to zero psi, and cavitation occurs. Once cavitation 
occurs, the barge is released from the fluid, and the structural velocity increases at a 
lower magnitude. The process continues again with the cavitating fluid particles pulling 
downward on the structure and releasing the structure until this cycle eventually settles 
out. 

Damping effects are present within the fluid and structural nodes as a result from 
the initial excitation (explosion) into the fluid medium. In LS-DYNA’s 
*CONTROL_ALE command, the Rayleigh damping was set to the on position. The 
Gruneisen EOS takes into account damping effects within the fluid particles. The fluid 
particle’s damping effect is then imparted to the structural nodes and produces an overall 
smooth velocity response. Similarities exist with previous research conducted 
concerning Rayleigh damping effects through the use of the LS-DYNA/USA program 
code [Ref. 9]. 
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Node 

Numerical Result 

Pmax (pSl) 

Empirical Result 

Pmax (psi) 

% Error 

1128 

1682.0 

2559.7 

34.0 

1139 

2304.3 

3113.2 

26.0 

46247 

4998.4 

5284.5 

5.4 

46256 

3387.2 

3590.6 

5.7 

46264 

2369.0 

2894.1 

18.1 


Table 5. Rectangular Block Case Results: Nodes Along X-Y Plane 
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Node 

Numerical Result 

P max (psi) 

Empirical Result 

Pmax (psi) 

% Error 

67155 

2137.9 

2694.7 

21.0 

69178 

1793.5 

2711.7 

33.9 

79616 

1790.1 

2648.5 

33.8 

79786 

1495.9 

2416.9 

38.1 

89068 

1105.1 

2022.24 

45.4 


Table 6. Rectangular Block Case Results: Nodes Along NRB 1 
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Node 

Numerical Result 

P max (psi) 

Empirical Result 

Pmax (psi) 

% Error 

25576 

1086.4 

1930.6 

43.7 

43680 

1394.8 

2150.1 

35.1 

46087 

2395.8 

2453.7 

2.4 

47694 

1481.5 

2222.1 

33.3 

47706 

1106.5 

1956.8 

43.4 


Table 7: Rectangular Block Case Results: Nodes Along NRB 2 
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Figure 35. Barge Nodes in Rectangular Box Charge Case 
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Figure 36. Barge Nodes in Stair Charge Case 
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Figure 37. Fluid Nodes in the Rectangular Box Charge Case 
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Figure 38. Fluid Nodes in the Stair Charge Case 
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Fluid Element 1128 Pressure Profile 



Figure 39. Fluid Element 1128 Pressure Profile 
(Rectangular Block Case Along XY Plane) 
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Figure 42. Fluid Element 46256Pressure Profile 
(Rectangular Block Case Along XY Plane) 
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Fluid Element 67155 Pressure Profile 
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Figure 44. Fluid Element 67155 Pressure Profile 
(Rectangular Block Case Along NRB 1) 




































Fluid Element 79616 Pressure Profile 



Figure 46. Fluid Element 79616 Pressure Profile 
(Rectangular Block Case Along NRB 1) 
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Fluid Element 79786 Pressure Profile 
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Figure 47. Fluid Element 79786 Pressure Profile 
(Rectangular Block Case Along NRB 1) 
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Figure 48. Fluid Element 89068 Pressure Profile 
(Rectangular Block Case Along NRB 1) 
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Figure 49. Fluid Element 25576 Pressure Profile 
(Rectangular Block Case Along NRB 2) 


























Fluid Element 43680 Pressure Profile 



Figure 50. Fluid Element 43680 Pressure Profile 
(Rectangular Block Case Along NRB 2) 
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Figure 55. Velocity Response at Nodes 124024/16862 
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Node 124173 vertical Velocity 



Figure 57. Velocity Response at Nodes 124173/39868 
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Node 124190 vertical veloci 



Figure 58. Velocity Response at Nodes 124190/40492 

















2 . 


Stair Explosive Simulation 


In the stair explosive simulation, sixteen solid fluid elements were used for 

comparison. Along the X-Y plane, nodes 3277, 3290, 73350, 73360, and 73370 were 
examined. Along NRB 1 plane, nodes 105852, 112932, 126792, 127012, and 141092 
were compared to the empirical formulations. And along the NRB 2 plane, nodes 40981, 
40981, 75223, 75245, and 75261 were used for comparison. The last fluid element 
involved 74296 was located underneath the barge’s centerline approximately. Figures 59 
through 74 illustrate the agreeable numerical and empirical pressure results. Tables 8 
through 10 compare the results in maximum peak pressure, P max , between the empirical 
and numerical formulations. The maximum pressures in this charge model are more 
agreeable than in the previous case. This is due to modeling the charge as close to a 
sphere than in the previous case of the rectangular box. The shock pressure appears to 
travel more like a spherical wave because of the stair charge geometry. 

Cavitation can also be observed in several nodes, namely nodes 73370, 75261, 
and 74296, which were at depths of 2.3 ft, 2.5 ft, and 1.0 ft respectively. The mentioned 
nodes at their respective depths were well within the cavitation zone and were in 
agreement with Figure 28. The other nodes were outside of the cavitational region. As 
expected, the water particle velocities at the fluid-air boundary behaved principally in a 
vertical manner due to the bulk cavitation. Additionally, the air pressure near the vicinity 
of the fluid-air boundary experienced little or no relative change in pressure as predicted. 

The structural nodes 195304, 195411, 195450, and 195467 correspond to the 
following fluid nodes, respectively, 26928, 61952, 61950, and 62654. In some instances, 
the structural node location did not match up exactly with the fluid node location, so the 
closest fluid node was chosen. Figures 75 through 78 illustrate the rapidly increasing 
vertical velocities (approximately 5 to 6 msec) to a peak value and then a rapid decrease 
and increase in a series of steps until the response settles down. The response does not 
settle out at a value of zero due to rigid body motion of the structure. Damping effects 
are also present as in the previous case. 
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In the rectangular box charge model, a total of 1 hour, 45 minutes was required 
for the simulation to complete its numerical computations. A time scale factor of 0.67 
was used in both simulation runs. The average time step was 12.5 (isec for 10 msec of 
data. In the stair charge model, a total of 2 hours, 50 minutes was required to run the 
simulation. The average time step was 12.5 jisec for 10 msec of data. Although LS- 
DYNA has the capability to run parallel processing, the option was not available for this 
numerical processing technique during the time of the experimentation. Once the parallel 
processing is fully incorporated in this ALE technique, the processing time will be 
significantly reduced. 

Overall, a good agreement in pressure signature shows the ALE technique is 
capable of handling the underwater explosion problem. A confidence in the numerical 
procedure has been established for this underwater explosion problem. A shift in the 
emphasis will now be focused on the two-dimensional air explosion problem interacting 
with a structural plate. 
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Node 

Numerical Result 

Pmax (psi) 

Empirical Result 

Pmax (psi) 

% Error 

3277 

2147.1 

2595.5 

17.3 

3290 

2964.9 

3240.4 

8.5 

73350 

6417.0 

5280.1 

21.5 

73360 

3798.2 

3740.5 

1.5 

73370 

2588.4 

2912.0 

11.1 


Table 8. Stair Case Results: Nodes Along X-Y Plane 




















Node 

Numerical Result 

Pmax (psi) 

Empirical Result 

Pmax (psi) 

% Error 

105852 

2565.0 

2693.3 

4.8 

112932 

2277.4 

2705.7 

15.8 

126792 

2204.1 

2643.8 

16.6 

127012 

1790.1 

2445.9 

26.8 

141092 

1261.2 

2012.5 

37.3 


Table 9. Stair Case Results: Nodes Along NRB 1 
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Node 


Numerical Result 

Pmax (PSI) 


Empirical Result 

Pmax (psi) 


% Error 






40981 

1380.2 

1937.7 

28.8 

40993 

1700.3 

2135.1 

20.4 

75223 

2797.5 

2453.7 

14.0 

75245 

1814.2 

2231.7 

18.8 

75261 

1275.7 

1948.7 

34.5 


Table 10: Stair Case Results: No4es Along NRB 2 

























Figure 59. Fluid Element 3277 Pressure Profile 
(Stair Case Along XY Plane) 
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Figure 60. Fluid Element 3290 Pressure Profile 
(Stair Case Along XY Plane) 
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Figure 61. Fluid Element 73350 Pressure Profile 
(Stair Case Along XY Plane) 
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Fluid Element 73360 Pressure Profile 
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Figure 62. Fluid Element 73360 Pressure Profile 
(Stair Case Along XY Plane) 
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Figure 63. Fluid Element 73370 Pressure Profile 
(Stair Case Along XY Plane) 
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Figure 64. Fluid Element 105852 Pressure Profile 
(Stair Case Along NRB 1) 
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Fluid Element 112932 Pressure Profile 
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Figure 66. Fluid Element 126792 Pressure Profile 
(Stair Case Along NRB 1) 
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Figure 67. Fluid Element 127012 Pressure Profile 
(Stair Case Along NRB 1) 
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Fluid Element 40981 Pressure Profile 



Figure 69. Fluid Element 40981 Pressure Profile 


(Stair Case Along NRB 2) 























Fluid Element 40993 Pressure Profile 



Figure 70. Fluid Element 40993 Pressure Profile 
(Stair Case Along NRB 2) 
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Figure 72. Fluid Element 75245 Pressure Profile 
(Stair Case Along NRB 2) 















Fluid Element 75261 Pressure Profile 



Figure 73. Fluid Element 75261 Pressure Profile 
(Stair Case Along NRB 2) 
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Fluid Element 74296 Pressure Profile 
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Figure 74. Fluid Element 74296 Pressure Profile 
(Stair Case Underneath Barge) 
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Node 195467 Vertical Velocity 



Figure 75. Velocity Response at Nodes 195467/62654 
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Node 195450 Vertical velocity 



Figure 76. Velocity Response at Nodes 195450/61950 
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Node 195411 vertical Velocity 



Figure 78. Velocity Response at Nodes 195411/61952 
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B. TWO-DIMENSIONAL AIR EXPLOSION 


As a result of the numerical agreement in the underwater explosion problem, the 
ALE formulation will be applied to a two-dimensional air explosion problem. It is the 
intent of this section to demonstrate and lay a foundation on the blast pressure wave 
signatures in the air medium as well as effects on a structural surface. Three simulations 
were observed for comparison by varying the structural thickness surrounding the TNT 
explosive. Although a 2.54 cm (1-inch) air gap existed between the explosive and 
structural cylinder, other distances could have been chosen. Furthermore, thermal effects 
were neglected due to the EOS was not capable of handling temperature effects. 

In Figures 79 through 81, the pressure rises quickly to maximum pressure, P max , at 
a certain time and returns to ambient pressure. The numerical pressure signatures 
simulated using the ALE technique show similar pressure responses found in earlier 
studies [Ref. 16]. Air pressure also decreases rapidly when the air element is farther 
away from the blast center. Notably, the structural cylinder’s thickness has a definite 
effect on the surrounding air pressure. As the thickness of the air cylinder increases, the 
maximum air pressure also decreases. 

The air particle velocities also exhibited a quick rise to a maximum velocity and 
experienced exponentially decay thereafter (Figures 82 through 84). The air particle 
velocities show that they reach their maximum velocities at the approximately the same 
time the maximum air pressure wave arrives. Furthermore, the air particle velocities also 
decrease at a rapid rate when the air particles are farther away from the blast center. 

As the explosion initiates, the cylindrical wall experiences deformation and 
ultimately ruptures. Although the material’s ultimate strength was not used, the structural 
shell elements were eroded through the material’s yield strength in LS-DYNA’s material 
card for the steel in order to provide the needed illustration effects for these models. 
More research is yet required is to find a suitable tangent modulus and hardening 
parameter as these factors play a role on material deformation and failure under an 
applied pressure. As noted in Figure 85, the structural velocities were notably high, 
especially for the thinnest shell elements. 
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Air Element 2515 Pressure Profile 



Figure 79. Air Pressure Signatures at Radius = 2.09 ft 
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Air Element 2528 Pressure Profile 



Figure 80. Air Pressure Signatures at Radius = 3.60 ft 
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Air Element 2555 Pressure Profile 



Figure 81. Air Pressure Signatures at Radius = 6.73 ft 
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Air Node 2292 Velocity Profile 



Figure 82. Air Velocity at Radius = 2.09 ft 
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Air Node 3904 velocity Profile 



Figure 83. Air Velocity at Radius = 3.60 ft 
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Shell Element 11043 Velocity Profile 
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VII. CONCLUSIONS AND RECOMMENDATIONS 


The Arbitrary Lagrangian-Eulerian (ALE) numerical technique was used to model 
a three-dimensional underwater explosion and two-dimensional air explosion. The 
simulations involved modeling the formation and propagation of the water and air shock 
waves as well as the fluid-structure and air-structure effects. In the underwater explosion 
shock wave propagation, the ALE’s numerical results were verified with empirical 
results. With satisfactory results, a numerical confidence was established in the 
underwater explosion problem. As a result, a two-dimensional air explosion was 
modeled. The air explosion simulation achieved similar characteristic results from 
previous studies, although a different processing program code was used in those 
simulations [Ref. 16]. 

Through material modeling and using various equations of state, the study 
incorporated numerous modeling facets. It is recommended that additional studies be 
conducted to examine further the extent of the ALE formulations. Specifically, the 
following areas require additional study: 

1. The effect of varying the penalty-coupling factor as it effects the fluid/air- 
structure interface. 

2. The effect of reducing the time scale factor i.e. achieving better numerical 
results versus computational time. 

3. Compare other different shock Hugoniots for water using the ALE numerical 
formulation. 

4. Investigate the strain effects in the plastic kinematic/isotropic material model. 

5. Investigate an equation of state, which takes into account the thermal effects 
for an air explosion. 
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APPENDIX A. BULK CAVITATION PROGRAM 


The following program code determines the bulk cavitation zone by solving 
Equations (3.18) and (3.19). The code is written for MATLAB Version 5.2 [Ref. 23]. 

% ME 4525 DESIGN PROJECT 1 
% by LT Ted Trevino 

% This program computes the axisymmetric upper and lower cavitation 
%boundaries for an 

% input charge (TNT or HBX-1). The program iterates up to a cavitation 
%zone radius 

% of 1400 ft and vertical distance of -50 ft. The program also 
%calculates the maximum 

% horizontal radius at its certain depth, and area enclosed by the upper 
%and lower 

% boundaries. Cavitation boundaries are based on Aron's Method. 

clear 
cl f 

gamma = .03703; % seawater weight density @ 1 atm (lbm/in A 3) 

pa = 14.7; % atmospheric pressure (psi) 

c = 5000; % acoustic velocity (ft/s) 

disp(' ') 

disp('Select charge type number: '); 
disp(' (1) HBX-1') 

disp ( ' (2) TNT') 

charge_type = input('Enter selection number: 

W = input('Enter charge weight (lbf): 

D = input('Enter charge depth (ft): 

% HBX-1 Calculations 
if (charge_type == 1) 

chg_name = 'HBX-1 Charger- 

weight - num2str(W); 

depth = num2str(D); 

K1 = 22347.6; 

K2 = 0.056; 

A1 = 1.144; 

A2 = -0.247; 

ub_data = [ ] ; 

lb_data = [ ] ; 

% Calculate upper boundary 
for x = 0:1400 

for y = 0:0.1:50; 

rl = sqrt((D-y) A 2 + x A 2) ; % Standoff distance charge-point 

r2 = sqrt((D+y) A 2 + x A 2); % Standoff distance image-point 

theta = K2*W A (1/3) * (W A (1/3)/rl) A A2/1000; % decay constant (sec) 
Fterml = (Kl*(W A (1/3)/rl) A Al*exp(-(r2-rl)/(c*theta))); 

Fterm2 = (gamma*y*12) - (Kl* (W A (1/3 )/r2) A A1) ; 


% Constant 
% Constant 
% Constant 
% Constant 

% Upper boundary data for x and y values 
% Lower boundary data for x and y values 


') ; 
') ; 
'); 
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% Upper cavitation boundary 


F = Fterml + Fterm2 + pa; 
if F <= 0 

ub__data = [ub_data; F x - (y) ] ; 
break 

end 

end 

end 

% Calculate lower boundary 
for x = 0:(length(ub_data)-1) 
for y = 0:0.1:50 

rl = sqrt((D-y) A 2 + x A 2) ; % Standoff distance charge-point 

r2 = sqrt((D+y) A 2 + x A 2); % Standoff distance image-point 

theta = K2*W A (1/3)*(W A (1/3)/rl) A A2/1000; % Decay constant 
% Incident SW pressure 

pi = Kl*(W A (1/3)/rl) A Al*exp(-(r2-rl)/(c*theta)); 

Gterml = - (pi / (c*theta)-) *(1+(((r2- 

(2*D*(D+y)/r2))/rl)*(((A2*r2)/rl)-A2-1))); 

Gterm2 = -((Al*pi)/rl A 2)*(r2-2*D*((D+y)/r2)); 

Gterm3 = (gamma*12)*((D+y)/r2); 

Gterm4 = (Al/r2)*(pi+pa+(gamma*y*12)); 

% Lower cavitation boundary 
G = Gterml + Gterm2 + Gterm3 + Gterm4; 
if G >= 0 

lb_data = [lb_data; G x -(y)]; 
break 

end 

end 

end 

% TNT Calculations 
elseif (charge_type == 2) 
chg_name = 'TNT Charge'; 
weight = num2str(W); 
depth = num2str(D); 

Kl = 22505; 

K2 = 0.058; 

A1 =1.18; 

A2 = -0.185; 

ub_data = [ ] ; 

lb_data = [ ] ; 

% Calculate upper boundary 
for x = 0:1400 

for y = 0:0.1:50; 

rl = sqrt((D-y) A 2 + x A 2); % Standoff distance charge-point 

r2 = sqrt((D+y) A 2 + x A 2); % Standoff distance image-point 

theta = K2*W A (1/3)*(W A (1/3)/rl) A A2/1000; % decay constant (sec) 
Fterml = (Kl*(W A (1/3)/rl) A Al*exp(-(r2-rl)/(c*theta))); 

Fterm2 = (gamma*y*12) - (Kl*(W A (1/3)/r2) A A1); 

F = Fterml + Fterm2 + pa; % Upper cavitation boundary 
if F <= 0 

ub_data = [ub_data; F x -(y)]; 
break 

end 


% Constant 
% Constant 
% Constant 
% Constant 

% Upper boundary data for x,y values 
% Lower boundary data for x,y values 
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end 


end 

% Calculate lower boundary 
for x = 0:(length(ub_data)-1) 
for y = 0:0.1:50 

rl = sqrt ( (D-y) A 2 + x A 2); % Standoff distance charge-point 

r2 = sqrt((D+y) A 2 + x A 2); % Standoff distance image-point 

theta = K2*W A (1/3)*(W A (1/3)/rl) A A2/1000; % Decay constant 
% Incident SW pressure 

pi = Kl*(W A (1/3)/rl) A Al*exp(-(r2-rl)/(c*theta)); 

Gterml = -(pi/(c*theta))*(1+(((r2- 

(2*D*(D+y)/r2))/rl)*(((A2*r2)/rl)-A2-1))) ; 

Gterm2 = -((Al*pi)/rl A 2)*(r2-2*D*((D+y)/r2)); 

Gterm3 = (gamma*12)*((D+y)/r2); 

Gterm4 = (Al/r2)*(pi+pa+(gamma*y*12)); 

% Lower cavitation boundary 
G = Gterml + Gterm2 + Gterm3 + Gterm4; 
if G >= 0 

lb__data = [lb_data; G x -(y)]; 
break 

end 

end 

end 
else 

disp('You did not enter selection 1 or 2. Please rerun program with 
appropriate entry.'); 
break 

end 

% Define boundary cavitation envelope 
ub = f ind(ub__data ( : , 3 ) <lb__data (:, 3) ) ; 
lb = find(lb_data(:,3)>ub_data(:,3)); 
ub_data(ub,:) = []; 

lb__data (lb, :) = []; 

height = ub_data(:,3) 
area = sum(height) ; 

vertl = [ -D] ; 

horl - [ 0]; 

%disp(' ') 

%disp f'****************************************************************** 
*★**★*■*•*★★★*■*•' ^ 

%fprintf('The maximum horizontal radius is %1.2f ft at a depth of %1.2f 
ft.\n',hor,vert) 

%disp( / ') 

%fprintf('The 2-D axisymmetric bulk cavitation envelope area = %1.2f 

ft A 2\n',area) 

orient landscape 

% Plot upper and lower boundary data and max horizontal radius 

plot(ub_data(:,2),ub_data(:,3),'b',lb_data(:,2),lb_data(:,3),'g ; ,horl,ver 

tl,' m* ' ) ; 

grid 


- lb_data(:,3);% Heights bw up & low boundaries 
% Area enclosed by upper/lower boundaries 
% Charge location vertically 
% Charge location horizontally 
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title(['Bulk Cavitation Envelope for ',weight,' lb ',chg_name,' at 
',depth,' ft']); 

xlabel('Horizontal Range (ft)') 
ylabel('Vertical Depth (ft)') 

legend ('Cavitation Upper BoundaryCavitation Lower Boundary', 4) 
axis([0 800 -55 5]) 
text(25,-2-D, 'Charge') 
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APPENDIX B. USEFUL TRUEGRID COMMANDS 


TrueGrid [Refs. 7 and 18] is a powerful finite element-modeling tool whose purpose 
is to create multiple block-structured mesh. The program has many useful features to aid in 
visualizing and manipulating a finite element model. Different input and output formats are 
also supported, such as LS-DYNA keyword format. Basic familiarity with TrueGrid is 
assumed. TrueGrid consists of three phases: control, part and merge phases. In the control 
phase, material models are defined and initialized. In the part phase, the mesh is 
constructed through projecting, deleting, and refining the mesh. Boundary conditions are 
also initialized in this phase. In the merge phase, the different parts are assembled together 
by merging nodes. Output formats (deck) are constructed during this phase. This appendix 
will cover some of the important commands used for mesh generation. 

A block or cylinder command can be entered in the control or merge phases. The 
block command syntax consists of i_indices; jjndices; kjndices; x_coordinates; 
y.coordinate; ^coordinates. Once the block command has been entered, the i, j, k indices 
are found in the computational window while the corresponding x, y, z coordinates exist in 
the physical window. The cylinder command syntax consists of ijndices; jjndices; 
kjndices; r_coordinates; 6_coordinates; z_coordinates. 

Once in the parts phase i.e. a block or cylinder command has been issued, a material 
number should be established through the mate command. Boundary conditions were 
established using the nodal displacement and rotation constraints, b, and non-reflecting 
boundaries, nr, commands. Elements from a part can be deleted using a combination of the 
computational window and environmental window. Elements needing to be deleted can be 
highlighted in the computational window using the mouse, then the delete button can be 
selected in the environmental window. This should generate a new mesh in the physical 
and computational windows with the elements deleted. For explosives, the detonation 
point, detp, command specifies where the detonation point will be placed in the actual 
physical window. Lastly, to prevent the structural nodes from merging with the fluid nodes 
a sliding interface region, si, was placed. The fluid nodes were assigned as the master 
nodes, and the structural nodes were assigned as the slave nodes. 
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Since the models deal with multi-materials or multi-parts, the nodes must be merged 
to form some type connectivity. Once all of the parts for the model were entered, the nodes 
were merged using the stp command. The nodes that fall within a specified tolerance 
defined by the user will merge with the exception of the sliding interface regions. Point 
masses are added using the pm command. After the entire mesh has been generated, LS- 
DYNA output format was selected using the OUTPUT command and then the WRITE 
command with no arguments is selected. This will produce an LS-DYNA format output 
file name called "trugrdo." 
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APPENDIX C. TRUEGRID COMMANDS: UNDEX 3-D MODEL WITH BLOCK 

TNT EXPLOSIVE MODEL 


This Appendix provides the input program used for the underwater three- 
dimensional model using a block TNT explosive model approach. The structural material 
properties and thickness values were input in LS-DYNA. TrueGrid does offer the feature 
to input material properties using the LSDYMATS command. Appendix B and references 
7 and 28 offer the meaning to the command lines. 


c ted trevino 

c tnt block explosive model 

c establish explosive mesh using block command 

block 123 4;1 2;1 2;-27.405 -9.135 9.135 27.405;0 17.45;0 17.45; 
c set up axis-symmetric boundaries 
b 1 1 1 4 1 2 dy 1 ry 1 ; 
b 1 1 1 4 2 1 dz 1 rz 1; 

c explosive material number that will be used to describe nodes in 
c LS-DYNA. 

c Material card can be defined in TrueGrid but is defined in LS-DYNA 
mate 11 

c Establish detonation point for the explosive 

detp 11 point -9.135 0 0; 

merge 

c water mesh 

block 1 18 29 30 31 32 43 60;1 2 41 43;1 2 18 22 39;-518.16 -213.36 - 
27.405 

-9.135 9.135 27.405 213.36 518.16;0 17.45 701.04 731.52;0 17.45 304.8 

365.76 

670.56; 

c delete unwanted water mesh to make room for barge and explosive 
dei 3 6; 1 2; 1 2; 

dei 2 7; 3 4; 3 4; 

c water material nr that will be used in LS-DYNA 
mate 10 

c Establish axis-symmetric boundaries 
b 1 1 1 8 4 1 dz 1 rz 1; 

b 1 1 1 8 1 5 dy 1 ry 1; 

c Establish non-reflecting boundaries 
nr 1 1 1 1 4 5 

nr 1 1 5 8 4 5 

nr 8 1 1 8 4 5 

c setup interface along the water nodes next to the barge so those nodes 
c will not merge with the barge. ! Important! Defined water as the master 
c nodes 

sid 1 lsdsi pi ; ; 
si 2332441m; 
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si2337341m; 
si7337441m; 
si 2337431m ; 
si 2347441m; 
merge 

c air mesh 

block 1 18 29 30 31 32 43 60;1 3 13;1 2 18 22 39;-518.16 -213.36 -27.405 
-9.135 

9.135 27.405 213.36 518.16;701.04 731.52 914.4;0 17.45 304.8 365.76 
670.56; 

c delete unwanted air nodes 

dei 1 8; 12; 13; 

dei 1 8; 1 2; 4 5; 

dei 7 8; 1 2; 3 4; 

dei 1 2; 1 2; 3 4; 

c establish axis-symmetric boundary 

b 1 1 1 8 3 1 dz 1 rz 1; 

c establish non-reflective boundaries 

nr 1 1 1 1 3 5 

nr 8 1 1 8 3 5 

nr 1 1 5 8 3 5 

c air material number 

mate 12 

merge 

c barge mesh 

block -1 -6 -12 -17;-1 7;-l 4 -7;-213.36 -86.36 86.36 213.36;701.04 

762;3 04.8 

335.28 365.76; 

c barge material number 

mate 1 

c establish interface along the barge shells so nodes will not merge with 

c with adjacent water. Defined as slave nodes. 

si 1111231s ; 

si 1114131s; 

si 4114231s ; 

si 1114211s; 

si 1134231s; 

merge 

c merge node tolerance 
stp 0.1 

c add point masses along three nodes 
pm 124126 1.16338e+5; 
pm 124261 1.1633 8e+5; 
pm 124392 1.16338e+5; 
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APPENDIX D. TRUEGRID COMMANDS: UNDEX 3-D MODEL WITH STAIR 

TNT EXPLOSIVE MODEL 


This Appendix provides the input program used for the underwater three- 
dimensional model using a stair-cased TNT explosive model approach. The structural 
material properties and thickness values were input in LS-DYNA. TrueGrid does offer the 
feature to input material properties using the LSDYMATS command. Appendix B and 
references 7 and 27 offer the meaning to the command lines. 


c ted trevino 

c tnt stair-cased charge model 
c charge model 

block 1 2 3 4 5 6;1 2 3 4;1 2 3 4;-26.5 -15.91 -5.3 5.3 15.91 26.5;0 
10.60 

21.2 31.80;0 10.60 21.2 31.80; 

c delete numerous blocks to model explosive stair case model 

dei 12; 14; 3 4; 

dei 2 3; 14; 3 4; 

dei 12; 1 4; 2 3; 

dei 1 2; 2 4; 1 2; 

dei 2 3; 3 4; 1 2; 

dei 2 3; 3 4; 2 3; 

dei 2 3; 2 3; 2 3; 

dei 3 4; 2 4; 3 4; 

dei 3 4; 3 4; 2 3; 

dei 46; 1 4; 3 4; 

dei 4 5; 2 4; 2 3; 

dei 5 6; 1 4; 2 3; 

dei 4 6; 3 4; 1 2; 

dei 5 6; 2 3; 1 2; 

c charge material number 

mate 11 

c axis symmetric boundary 
b 1 1 1 6 4 1 dz 1 rz 1 ; 

b 1 1 1 6 1 4 dy 1 ry 1 ; 

c establish detonation point 
detp 11 point -5.30 0 0; 
merge 

c fluid model 

block 1 21 33 34 35 36 37 38 50 70;1 2 3 4 48 50;1 2 3 4 21 25 45;-518.16 
-213.36 -26.50 -15.90 -5.30 5.30 15.90 26.50 213.36 518.16;0 10.60 21.20 
31.80 

701.04 731.52;0 10.60 21.20 31.80 304.8 365.76 670.56; 
c make deletions for fluid to make room for the explosive and barge 
dei 38; 1 2; 1 2; 

dei 47; 23; 12; 
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dei 4 7; 1 2; 2 3; 

dei 5 6; 3 4; 1 2; 

dei 5 6; 12; 3 4; 
dei 5 6; 2 3; 2 3; 

dei 2 9; 5 6; 5 6; 

c water material nr 
mate 10 

c axis-symmetric boundaries 
b 1 1 1 10 6 1 dz 1 rz 1 ; 

b 1 1 1 10 1 7 dy 1 ry 1 ; 

c non-reflective boundaries 
nr 1 1 1 1 6 7 
nr 10 1 1 10 6 7 
nr 1 1 7 10 6 7 

c establish interface along the water nodes that will interface on the 

c barge nodes. Master nodes. 

sid 1 lsdsi pi ; ; 

si 2552661m; 

si 2559561m ; 

si 9559661m; 

si2559651m; 

si2569661m; 

merge 

c air model 

block 1 21 33 34 35 36 37 38 50 70;1 3 15;1 2 3 4 21 25 45;-518.16 - 

213.36 

-26.50 

-15.90 -5.30 5.30 15.90 26.50 213.36 518.16;701.04 731.52 914.4;0 10.60 
21.20 

31.80 304.80 365.76 670.56; 
c delete air model portions 
dei 1 10; 1 2; 1 5; 
dei 1 10; 12; 67; 
dei 1 2; 12; 56; 
dei 9 10; 12; 56; 
c air material number 
mate 12 

c axis symmetric boundaries 

b 1 1 1 10 3 1 dz 1 rz 1; 

c nonreflective boundaries 

nr 1 1 1 1 3 7 

nr 10 1 1 10 3 7 

nr 1 1 7 10 3 7 

merge 

c structure model 

block -1 -6 -12 -17;-1 7;-l 4 -7;-213.36 -86.36 86.36 213.36;701.04 
762;304.8 
335.28 365.76; 
c barge material nr 
mate 1 

c interface on the structure. Slave nodes 
si 1111231s ; 
si 1114131s; 
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si4114231s; 
si 1114211s ; 
si 1134231s; 
merge 

c merge tolerance 
stp 0.1 

c point masses 
pm 195338 1.16338e+5 
pm 195471 1.16338e+5 
pm 195607 1.16338e+5 






THIS PAGE INTENTIONALLY LEFT BLANK 
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APPENDIX E. LS-DYNA INPUT DECK: UNDEX 3-D MODEL WITH BLOCK TNT 

EXPLOSIVE MODEL 


This Appendix provides the input deck used for the underwater three-dimensional 
model using a block TNT explosive model approach. Appendix G and references 13 and 
25 offer the meaning to the command lines. 


$ Barge model with TNT explosive offset. 

$ Ted Trevino 
*KEYWORD 
*TITLE 

Barge TNT Coupled Eulerian/Lagrangian 

* C ONTRO L_TERMINATION 

10000.,0,0.000E+00,O.OOOE+OO,O.OOOE+OO 

*CONTROL_TIME_STEP 

0,0.67 

* DATABAS E_BINARY_D3PLOT 

200 . 0,0 

*CONTROL_ENERGY 

2 , 1 , 2,2 

*DATABASE_BINARY_D3 THDT 

1 . 0,0 

*DATABASE_HISTORY_NODE 

124190,124173,124024,124139,124149,124134,124164,124287 

* DAT ABA S E_NODOUT 
1.0 

* DATABAS E_GLSTAT 
1. 

* CONTROL_ALE 

2 , 1 , 2,-1 

0. 

*CONSTRAINED_LAGRANGE_IN_SOLID 
1,2,0,0,0,4,1,0 
0 , 0 , 0.1 

* SET_PART_LIST 
1 

1 

* SET_PART_LIST 
2 

10,12 

* ALE_MULTI -MATERI AL_GROUP_PART 
10 

* ALE_MULTI -MATERI AL_GROUP_PART 
11 

* ALE_MULTI -MATERI AL_GROUP_P ART 
12 

*ALE_INTERFACE_RECONSTRUCTION 

1 

$ 
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$ DEFINITION OF MATERIAL 1 (SHIP) 

$ 

*MAT_PLASTIC_KINEMATIC 

1,7.87,2.1,0.300,4.300E-03,0.000E+00,0.OOOE+OO 

0.000E+00,0.OOOE+OO,0.000E+00 

♦HOURGLASS 

1,0,0.OOOE+OO,0,0.OOOE+OO,0.OOOE+OO 
*SECTION_SHELL 

I, 2,0.OOOE+OO,3.00,0.OOOE+OO,0.OOOE+OO, 0 
0.6350,0.6350,0.6350,0.635,0.OOOE+OO 
♦PART 

material type # 3 (Elastic) 

1 , 1 , 1 , 0 , 0,0 

$ 

$ DEFINITION OF MATERIAL 10 (WATER) 

$ 

*MAT_NULL 

10,1.00,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 

* EO S_GRUNEISEN 

10,0.148,1.75,0.OOOE+OO,0.OOOE+OO,0.4934,0.OOOE+OO, 0.OOOE+OO 

1.00 

* SECTION_SOLID_ALE 

10,11,0 

0 

♦PART 

material type # 9 (Fluid) 

10 , 10 , 10 , 10 , 11,0 

$ 

$ DEFINITION OF MATERIAL 11 (TNT) 

$ 

*MAT_HIGH_EXPLOSIVE_BURN 

II, 1.63,0.693,0.210,0.OOOE+OO 
♦HOURGLASS 

11, 0,0,0,0,0,0 
*EOS_JWL 

11,3.71,3.230E-02,4.15,0.950,0.300,4.300E-02,1.0 

* SECTION_SOLID_ALE 

11 , 11,0 

0 

♦PART 

material type # 8 (High Explosive Burn) 

11 , 11 , 11 , 11 , 11,0 

$ 

$ DEFINITION OF MATERIAL 12 (AIR) 

$ 

*MAT_NULL 

12,1.280E-03,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 

* EOS_LINEAR_POLYNOMIAL 

12,0.00OE-O0,O.O,0.000E+00,O.OOOE+OO,0.400,0.400,0.OOOE+OO 
0.OOOE+OO,0.OOOE+OO 
*SECTION_SOLID_ALE 
12,11,0 
0 

♦PART 

material type # 9 (Fluid) 
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12,12,12,12,11,0 

$ 

$ NODES 

$ 

*NODE 

1, -27.4050007,0.OOOOOOOOOE+OO,0.000000000E+00,5,5 

2, -27.4050007,0.000000000E+00,17.4500008,2,2 

3, -27.4050007,17.4500008,0.OOOOOOOOOE+OO, 3,3 


124461,213.360001,751.840027,355.600037,0,0 
124462,213.360001,762.000000,345.440002,0,0 
124463,213.360001,762.000000,355.600006,0,0 
$ 

$ ELEMENT CARDS FOR SOLID ELEMENTS 

$ 

*ELEMENT_SOLID 

1,11,1,5,7,3,2,6,8,4 

2,11,5,9,11,7,6,10,12,8 

3,11,9,13,15,11,10,14,16,12 


116582,12,123492,123662,123679,123509,123493,123663,123680,123510 

116583,12,123662,123832,123849,123679,123663,123833,123850,123680 

116584,12,123832,124002,124019,123849,123833,124003,124020,123850 

$ 

$ ELEMENT CARDS FOR SHELL ELEMENTS 

$ 

*ELEMENT_SHELL_THICKNESS 

1,1,124021,124025,124026,124022 

0.635,0.635,0.635,0.635 

2,1,124025,124029,124030,124026 

0.635,0.635,0.635,0.635 

3,1,124029,124033,124034,124030 

0.635,0.635,0.635,0.635 


431.1.124459.124461.124432.124431 
0.635,0.635,0.635,0.635 

432.1.124461.124463.124433.124432 
0.635,0.635,0.635,0.635 

$ 

$ SLIDING INTERFACE DEFINITIONS 
$ 

$ 

$ TrueGrid Sliding Interface # 1 

$ 

* SET_SEGMENT 

1,0.OOOE+OO,0.000E+00,0.000E+00,0.000E+00 

124021,124025,124026,124022,0.000E+00,0.000E+00,0.000E+00,0.000E+00 
124022,124026,124027,124023,0.OOOE+OO,0.000E+00,0.000E+00,0.000E+00 
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124448,124451,124462,124460,0.000E+00,0.000E+00,0.OOOE+OO,0.OOOE+OO 
124460,124462,124463,124461,0.000E+00,0.000E+00,0.000E+00,0.000E+00 
124461,124463,124433,124432,0.000E+00,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 
* SET_SEGMENT 

2,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 

14056,29430,29583,16861,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 
16861,29583,29584,16862,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 
16862,29584,29585,16863,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 


99757,99760,99761,99758,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 
99758,99761,99762,99759,0.OOOE+OO, 0.000E+00,0.OOOE+OO, 0.OOOE+OO 
99759,99762,71944,71926,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 
$ 

$ DETONATION POINTS 
$ 

*INITIAL_DETONATION 

11,-9.14,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 
$ 

$ DISCRETE LUMPED MASSES 
$ 

*ELEMENT_MASS 

1,124126,116338.0 

2,124261,116338.0 

3,124392,116338.0 

$ 

$ NONREFLECTING BOUNDARY SEGMENTS 
$ 

* SET_SEGMENT 
3,0.,0.,0. ,0. 

17,19,20,18,0.,0.,0. ,0. 

18,20,105,89,0.,0. ,0. ,0. 


123799,123969,123986,123816,0.,0. ,0. ,0. 
123816,123986,124003,123833,0.,0.,0.,0. 
123833,124003,124020,123850,0.,0.,0.,0. 
* BOUNDARY_NON_REFLECTING 
3,0.,0. 

*END 
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APPENDIX F. LS-DYNA INPUT DECK: UNDEX 3-D MODEL WITH STAIR TNT 

EXPLOSIVE MODEL 


This Appendix provides the input deck used for the underwater three-dimensional 
model using a stair-cased TNT explosive model approach. Appendix G and references 13 
and 25 offer the meaning to the command lines. 


*KEYWORD 
*TITLE 

Barge TNT Coupled Eulerian/Lagrangian 
*CONTROL_TERMINATION 

100000,0.000E+00,0.000E+00,0.000E+00 

* CONTROL_TIME_STEP 
0,0.67 

*DATABASE_BINARY_D3 PLOT 

200 . 0,0 

* CONTROL_ENERGY 

2 , 1 , 2,2 

*DATABASE_BINARY_D3THDT 

1 . 0,0 

* DAT ABAS E_H I STORY_NODE 

195467,195304,195450,195416,195426,195411,195561,195441 

* DATABASE_NODOUT 
1.0 

*DATABASE_GLSTAT 

1 . 

* CONTROL_ALE 

2 , 1 , 2,-1 

0. 

*CONSTRAINED_LAGRANGE_IN_SOLID 
1,2,0,0,0,4,1,0 
0 , 0 , 0.1 

* SET_PART_LIST 
1 

1 

* SET_PART_LIST 
2 

10,12 

*ALE_MULTI-MATERIAL_GROUP_PART 

10 

*ALE_MULTI-MATERIAL_GROUP_PART 

11 

*ALE_MULTI-MATERIAL_GROUP_PART 

12 

*ALE_INTERFACE_RECONSTRUCTION 

1 

$ 

$ DEFINITION OF MATERIAL 1 (SHIP) 

$ 
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*MAT_PLASTIC_KINEMATIC 

1,7.87,2.1,0.300,4.300E-03,0.000E+00,0.000E+00 
0.000E+00,0.OOOE+OO,0.000E+00 
*HOURGLASS 

1,0,0.OOOE+OO,0,0.000E+00,0.OOOE+OO 
*SECTION_SHELL 

I, 2,0.OOOE+OO,3.00,0.OOOE+OO,0.OOOE+OO, 0 
0.6350,0.6350,0.6350,0.635,0.OOOE+OO 
♦PART 

material type # 3 (Elastic) 

1 , 1 , 1 , 0 , 0,0 

$ 

$ DEFINITION OF MATERIAL 10 (WATER) 

$ 

*MAT_NULL 

10,1.00,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 

* EOS_GRUNEISEN 

10,0.148,1.75,0.OOOE+OO,0.OOOE+OO,0.4934,0.OOOE+OO,0.OOOE+OO 

1.00 

*SECTION_SOLID_ALE 

10,11,0 

0 

♦PART 

material type # 9 (Fluid) 

10 , 10 , 10 , 10 , 11,0 

$ 

$ DEFINITION OF MATERIAL 11 (TNT) 

$ 

* MAT_HIGH_EXPLOSIVE_BURN 

II, 1.63,0.693,0.210,0.OOOE+OO 
♦HOURGLASS 
11 , 0 , 0 , 0 , 0 , 0,0 

*EOS_JWL 

11,3.71,3.230E-02,4.15,0.950,0.300,4.300E-02,1.0 
*SECTION_SOLID_ALE 
11 , 11,0 
0 

♦PART 

material type # 8 (High Explosive Burn) 

11 , 11 , 11 , 11 , 11,0 

$ 

$ DEFINITION OF MATERIAL 12 (AIR) 

$ 

*MAT_NULL 

12,1.280E-03,0.000E+00,0.OOOE+OO,0.000E+00,0.OOOE+OO 

* EOS_LINEAR_POLYNOMIAL 

12,0.OOOE-OO,0.0,0.OOOE+OO,0.OOOE+OO,0.400,0.400,0.OOOE+OO 
0.OOOE+OO,0.OOOE+OO 
*SECTION_SOLID_ALE 
12,11,0 
0 

♦PART 

material type # 9 (Fluid) 

12,12,12,12,11,0 

$ 
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$ NODES 

$ 

*NODE 

1, -26.5000000,O.OOOOOOOOOE+00,O.OOOOOOOOOE+O0,5,5 

2, -26.5000000,0.OOOOOOOOOE+OO,10.6000004,2,2 

3, -26.5000000,10.6000004,0.000000000E+00,3,3 


195731,213.360001,751.840027,355.600037,0,0 
195732,213.360001,762.000000,345.440002,0,0 
195733,213.360001,762.000000,355.600006,0,0 
$ 

$ ELEMENT CARDS FOR SOLID ELEMENTS 

$ 

* ELEMENT_SOLID 
1,11,1,5,7,3,2,6,8,4 
2,11,5,9,11,7,6,10,12,8 
3,11,6,10,12,8,13,15,16,14 


185194,12,194559,194799,194819,194579,194560,194800,194820,194580 
185195,12,194799,195039,195059,194819,194800,195040,195060,194820 
185196,12,195039,195279,195299,195059,195040,195280,195300,195060 
$ 

$ ELEMENT CARDS FOR SHELL ELEMENTS 

$ 

* ELEMENT_SHELL_THICKNES S 

1,1,195301,195305,195306,195302 

0.635,0.635,0.635,0.635 

2,1,195305,195309,195310,195306 

0.635,0.635,0.635,0.635 

3,1,195309,195313,195314,195310 

0.635,0.635,0.635,0.635 


430.1.195727.195729.195703.195702 
0.635,0.635,0.635,0.635 

431.1.195729.195731.195704.195703 
0.635,0.635,0.635,0.635 

432.1.195731.195733.181580.195704 
0.635,0.635,0.635,0.635 

$ 

$ SLIDING INTERFACE DEFINITIONS 
$ 

$ 

$ TrueGrid Sliding Interface # 1 

$ 

* SET_SEGMENT 

1,0.000E+00,0.000E+00,O.O00E+00,0.000E+00 

195301,195305,195306,195302,0.000E+00,0.000E+00,0.000E+00,0.OOOE+OO 
195302,195306,195307,195303,0.000E+00,0.000E+00,0.000E+00,0.OOOE+OO 
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195719,181578,195732,195730,O.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 
195730,195732,195733,195731,0.OOOE+OO,0.OOOE+OO,0.000E+00,O.000E+O0 
195731,195733,181580,195704, O.OOOE+OO, 0.000E+00,0.000E+00,0.000E+00 
* SET_SEGMENT 

2, O.OOOE+OO,O.OOOE+OO,O.OOOE+OO,O.OOOE+OO 

23234,46275,46453,26927,O.OOOE+OO,O.OOOE+OO,O.OOOE+OO,O.OOOE+OO 
26927,46453,46454,26928,O.OOOE+OO,O.OOOE+OO,O.OOOE+OO,O.OOOE+OO 


156367,156370,156371,156368,O.OOOE+OO,O.OOOE+OO, O.OOOE+OO,O.OOOE+OO 
156368,156371,156372,156369,O.OOOE+OO, O.OOOE+OO,O.OOOE+OO,O.OOOE+OO 
156369,156372,112306,112285,O.OOOE+OO,O.OOOE+OO, O.OOOE+OO,O.OOOE+OO 
$ 

$ DETONATION POINTS 

$ 

*INITIAL_DETONATION 

11,-5.30,O.OOOE+OO,O.OOOE+OO,O.OOOE+OO 

$ 

$ DISCRETE LUMPED MASSES 

$ 

*ELEMENT_MASS 

1,195338,116338.0 

2,195471,116338.0 

3,195607,116338.0 

$ 

$ NONREFLECTING BOUNDARY SEGMENTS 

$ 

* SET_SEGMENT 
3,0.,0.,0.,0. 

51,53,54,52,0.,0.,0.,0. 

52,54,136,135,0. ,0.,0.,0. 


195020,195260,195280,195040,0.,0.,0.,0. 
195040,195280,195300,195060,0.,0.,0.,0. 
* BOUNDARY_NON_REFLECTING 
3,0.,0. 

*END 
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APPENDIX G. USEFUL LS-DYNA COMMANDS 


LS-DYNA [Ref. 25] version 960 (alpha) was used for the numerical computations. 
This Appendix covers some of the basic commands used in the input decks. The pre- 
processor.(TrueGrid) primarily outputs into the "trgrdo" file a listing of nodes, solids, 
shells, boundary conditions, point masses, and detonation points. Note: material type and 
equation of state information could have been input in TrueGrid, but the data input was not 
in this case. Since the initial material information and equation of state was not done in the 
pre-processor, the data was input in LS-DYNA. Reference 25 provides a detailed 
description of the "keyword" commands. Some of the major keywords are listed below: 

a. *CONTROL_ALE sets default control parameters for the Arbitrary 
Lagrange-Eulerian and Eulerian calculations. This command works in 
conjunction with *ALE_MULTI-MATERIAL_GROUP, and 
*SECTION_SOLID_ALE. 

b. *ALE_MULTI-MATERIAL_GROUP_PART defines the PART ID s of 
each group pertaining to the Arbitrary Lagrangian-Eulerian capability. 

c. *SECTION_SOLID_ALE defines section material properties for solid 
structural and fluid elements. 

d. *CONSTRAINED_LAGRANGE_IN_SOLID couples the Lagrangian 
mesh (slaves) which can consist of shells, solids, or beams to an Eulerian 
flow (masters). In this command, the penalty coupling is initialized and set 
to a default of 0.10, although can be changed. May be subject to further 
study between in structural-fluid interface. 

e. *CONTROL_ENERGY provides controls for energy dissipation options. 

In this command, Rayleigh damping is initialized. 

f. *CONTROL_TIME_STEP sets the time step size control using a time 
scale factor of 0.67 or less for explosives. Any value higher may cause 
instabilities in the numerical solutions. Reducing the time scale factor may 
also may subject to additional studies. 
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g. *DATABASE_option controls database definitions that are necessary in 
obtaining output files containing results. Three-dimensional plotting 
information is contained which is used in LS-POST. ASCII information can 
also be obtained for nodal or element data which can also be used in LS- 
POST. 
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APPENDIX H. USEFUL LS-POST COMMANDS 


This Appendix contains some of the important LS-POST commands [Ref. 26]. LS- 
POST is a powerful three-dimensional post-processor, which can be used to display 
animation graphics as well element or nodal XY plots. In order to access the geometrical 
two or three-dimensional data, the d3plot must be selected. 

In order to view the shock pressure waves, the Fcomp button was selected in the 
Main Menu Area. A list of options will appear in the Fringe component window. The 
Misc button was then selected followed by the Pressure option. The Apply button is then 
selected in the Fringe component window. Lastly, to view to simulation run time pressure 

animation, the animation control buttons are used normally in the play position, ►. The 
other animation control buttons can manipulated to view the animation in reverse, one state 
at a time, and the animation speed can also be controlled. Most of the animation feature 
options or selections can be selected and applied from the Fcomp button such as stress, 
strain, and energy three-dimensional plots to name a few. 

To view the two or three-dimensional velocity vectors, the Vector button is selected 
from the Main menu area. Next, the velocity button is selected, and the vector length can 
also be scaled using the SF button to the user’s appropriate needs. Then, to run the velocity 
animation, follow the steps from the previous paragraph. 

To make a movie from the animation sequence that has been run in the graphics 
window, the user should select the File pull down menu and activate the Movie selection. 
The format used was AVI(rle) with the size on NTSC. The start button was then selected, 
and the animation sequence is recorded. 

XY plots can also be produced using the History button in the Main Menu Area. 
The user then selects either the Nodal or Element buttons in the Time History Results 
window. Selecting the Nodal button provides nodal coordinates, displacement, velocity 
and acceleration. Selecting the Element button provides the element stresses and strains. 
After the Nodal or Element selection has been made, the user then selects the actual node 
or element in the graphics window with the mouse. After the node or element has been 
chosen in the graphics window, the Plot button in the Time History Results window is 
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selected. The XY plot is then produced. The plot can be saved while in XY Graph window 
by selecting the Save button. This saves the XY data into file where the data is in a two- 
column format. 

Additionally, nodes or elements selected in LS-DYNA’s *database_option can be 
plotted using the Ascii button in the Main Menu Area. An Ascii File operation window 
will appear. The Load button was selected followed by the +nodout*. This generated a 
Nodout Data window from which the nodes that were selected from LS-DYNA’s 
*database_option will appear. The appropriate node was selected with for example Y-vel 
selected. The Plot button is then selected in the Nodout Data window. The XY plot is then 
produced. The plot can be saved using the steps from the previous paragraph. 

Lastly, the graphical image in the Graphics window of certain state can be saved 
into a file using the File pull down menu. Next, the Print command was selected. A Print 
Dialog window will appear and the print to File was selected. The file format selected was 
changed from the default of PS/Image to BMP, which produced the best results. 
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APPENDIX I. TRUEGRID COMMANDS: 2-D AIR EXPLOSION 


This Appendix provides the input program used for the two-dimensional air 
explosion using TNT. This program provides and generates the initial input mesh. All of 
the structural material properties and thickness values were input into LS-DYNA. These 
commands offer the template for the other two air explosion scenarios. TrueGrid does offer 
the feature to input material properties using the LSDYMATS command. Appendix B and 
references 7 and 25 offer the meaning to the command lines. 


c establish the air/explosive mesh 

cylinder 1 2 89;1 62;1 2;0 3.54 311.52;0 90;0 3.54; 

c air material nr. The nodes for the explosive are changed in LS-DYNA 
c to prevent using the interface command 
mate 12 

c setup the detonation point 
detp 11 point 000; 

c esabtlish axis-symmetric boundaries 
b 1 1 1 3 1 2 dy 1 ry 1 ; 
b 1 2 1 3 2 2 dx 1 rx 1 ; 

c constrain nodes in the z direction to maintain problem physics 

b 1 1 1 3 2 2 dz 1 rz 1 ; 

merge 

c establish cylindrical shell 
cylinder -1;1 30;1 2;6.08;0 90;0 3.54; 
c shell material number 
mate 1; 
merge 
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APPENDIX J. LS-DYNA INPUT DECK: 2-D AIR EXPLOSION WITH ONE-INCH 

THICK STEEL 


This Appendix provides the input deck used for the air explosion two-dimensional 
air explosion interacting with a one-inch thick steel. Appendix K and references 13 and 25 
offer the meaning to the command lines. 


*KEYWORD 
*TITLE 

APM Steel-TNT Coupled Eulerian/Lagrangian 

* CONTROL_TERMINATION 

4000.,0,0.000E+00,0.OOOE+OO,0.OOOE+OO 

*CONTROL_TIME_STEP 

0,0.67 

*DATABASE_BINARY_D3 PLOT 
50.0,0 

* CONTROL_ENERGY 

2 , 1 , 2,2 

*DATABASE_BINARY_D3 THDT 

1 . 0,0 

* DATABASE_HISTORY_NODE 
2000 

* DATABASE_NODOUT 

.8 

*DATABASE_GLSTAT 

1. 

*CONTROL_ALE 
2 , 1 , 2,-1 
0. 

* CONSTRAINED_LAGRANGE_IN_SOLID 
1,2,0,0,0,4,1,0 

0 , 0 , 0.1 

* SET_PART_L1ST 
1 

1 

* SET_PART_L1ST 

2 

10,12 

*ALE_MULTI-MATERIAL_GROUP_PART 

10 

* ALE_MULTI -MATERIAL_GROUP_PART 

11 

* ALE_MULTI -MATERI AL_GROUP_PART 
12 

*ALE_INTERFACE_RECONSTRUCTION 

1 

$ 
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$ DEFINITION OF MATERIAL 1 (Steel) 

$ 

*MAT_PLASTIC_KINEMATIC 

1,7.87,2.1,0.300,4.300E-03,0.OOOE+OO,0.000E+00 
0.OOOE+OO,0.000E+00,0.0 
*HOURGLASS 

1,0,0.OOOE+OO,0,0.000E+00,0.OOOE+OO 

* SECTION_SHELL 

I, 2,0.OOOE+OO,3.00,0.000E+00,0.000E+00,0 
2.54,2.54,2.54,2.54,0.OOOE+OO 

*PART 

material type # 3 (Kinematic/Isotropic Elastic-Plastic) 

1 , 1 , 1 , 0 , 0,0 

$ 

$ DEFINITION OF MATERIAL 10 (WATER) Not used in this simulation 

$ 

*MAT_NULL 

10,1.00,0.000E+00,0.000E+00,0.OOOE+OO,0.000E+00 

* EOS_GRUNEISEN 

10,0.148,1.75,0.000E+00,0.000E+00,0.280,0.000E+00,0.000E+00 

1.00 

*SECTION_SOLID_ALE 

10,11,0 

0 

*PART 

material type # 9 (Fluid) 

10 , 10 , 10 , 10 , 11,0 

$ 

$ DEFINITION OF MATERIAL 11 (TNT) 

$ 

*MAT_HIGH_EXPLOSIVE_BURN 

II, 1.63,0.693,0.210,0.000E+00 
*HOURGLASS 
11 , 0 , 0 , 0 , 0 , 0,0 

*EOS_JWL 

11,3.71,3.230E-02,4.15,0.950,0.300,4.300E-02,1.0 

* SECTION_SOLID_ALE 

11 , 11,0 

0 

*PART 

material type # 8 (High Explosive Burn) 

11 , 11 , 11 , 11 , 11,0 

$ 

$ DEFINITION OF MATERIAL 12 (AIR) 

$ 

*MAT_NULL 

12,1.280E-03,0.000E+00,0.OOOE+OO,0.OOOE+OO,0.000E+00 

* EO S_LINEAR_POLYNOMIAL 

12,0.OOOE-OO,0.0,0.000E+00,0.000E+00,0.400,0.400,0.000E+00 
0.000E+00,0.000E+00 

* SECTION_SOLID_ALE 

12,11,0 

0 

*PART 

material type # 9 (Fluid) 
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12,12,12,12,11,0 

$ 

$ NODES 

$ 

*NODE 

1,0.OOOOOOOOOE+OO,0.OOOOOOOOOE+OO,0.OOOOOOOOOE+OO,5,5 
2,0.OOOOOOOOOE+OO,0.OOOOOOOOOE+OO,3.53999996,5,5 
3,0.OOOOOOOOOE+OO,0.OOOOOOOOOE+OO,0.OOOOOOOOOE+OO,3,3 


11094,0.329162806,6.07108307,3.53999996,0,0 

11095,0.OOOOOOOOOE+OO,6.07999992,0.OOOOOOOOOE+OO,0,0 

11096,0.OOOOOOOOOE+OO,6.07999992,3.53999996,0,0 

$ 

$ ELEMENT CARDS FOR SOLID ELEMENTS 

$ 

*ELEMENT_SOLID 

1,11,1,125,127,3,2,126,128,4 

2,11,3,127,129,5,4,128,130,6 

3,11,5,129,131,7,6,130,132,8 


5367,12,10785,10909,10911,10787,10786,10910,10912,10788 
5368,12,10909,11033,11035,10911,10910,11034,11036,10912 
$ 

$ ELEMENT CARDS FOR SHELL ELEMENTS 

$ 

* EL EMENT_SHELL_THICKNES S 

1,1,11037,11039,11040,11038 

02.54,2.54,2.54,2.54 

2,1,11039,11041,11042,11040 

02.54,2.54,2.54,2.54 


29,1,11093,11095,11096,11094 

02.54,2.54,2.54,2.54 

$ 

$ DETONATION POINTS 

$ 

*INITIAL_DETONATION 

11,0.000E+00,0.000E+00,O.OOOE+OO,O.OOOE+OO 
♦END 
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APPENDIX K. LS-DYNA INPUT DECK: 2-D AIR EXPLOSION WITH HALF 

INCH THICK STEEL 


This Appendix provides the input deck used for the air explosion two-dimensional 
air explosion interacting with a half-inch thick steel. Appendix G and references 13 and 25 
offer the meaning to the command lines. 


*KEYWORD 
*TITLE 

APM Steel-TNT Coupled Eulerian/Lagrangian 

* CONTROL__TERMINAT I ON 

4000.,0,0.000E+00,0.000E+00,O.OOOE+OO 

*CONTROL__TIME_STEP 

0,0.67 

* DATAB AS E__B INARY_D 3 PLOT 
50.0,0 

*CONTROL_ENERGY 
2,1,2,2 

* DATABASESINARY_D3 THDT 

1 . 0,0 

* DATABA S E_HISTORY_NODE 
2000 

* DATABASE_NODOUT 

.8 

* DATABAS E_GL S TAT 
1. 

*CONTROL_ALE 
2 , 1 , 2,-1 
0. 

*CONSTRAINED_LAGRANGE_IN_SOLID 
1,2,0,0,0,4,1,0 
0 , 0 , 0.1 

* S ET_PART_L1ST 

1 

1 

* S ET__PART_L 1ST 
2 

10,12 

*ALE„MULTI-MATERIAL_GROUP__PART 

10 

*ALE_MULTI -MATERIAL__GROUP_PART 
11 

*ALE__MULTI-MATERIAL_GROUP_PART 

12 

*ALE_INTERFACE_RECONSTRUCTION 

1 

$ 
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$ DEFINITION OF MATERIAL 1 (Steel) 

$ 

*MAT__PLAST I C_KI NEMATIC 

1,7.87,2.1,0.300,4.300E-03,0.OOOE+OO,0.000E+00 
0.000E+00,0.OOOE+OO,0.0 
*HOURGLASS 

1,0,0.000E+00,0,0.OOOE+OO,0.000E+00 
*SECTION„SHELL 

1,2,0.OOOE+OO,3.00,0.OOOE+OO,0.OOOE+OO,0 

I. 27,1.27,1.27,1.27,0.OOOE+OO 
*PART 

material type # 3 (Kinematic/Isotropic Elastic-Plastic) 
1 , 1 , 1 , 0 , 0,0 
$ 

$ DEFINITION OF MATERIAL 10 (WATER) 

$ 

*MAT_NULL 

10,1.00,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO,0.OOOE+OO 
*EOS_GRUNEISEN 

10,0.148,1.75,0.OOOE+OO,0.OOOE+OO,0.280,0.OOOE+OO,0.OOOE+OO 

1.00 

* SECTION_SOLID_ALE 

10 , 11,0 

0 

*PART 

material type # 9 (Fluid) 

10 , 10 , 10 , 10 , 11,0 

$ 

$ DEFINITION OF MATERIAL 11 (TNT) 

$ 

*MAT_HIGH_EXPLOSIVE_BURN 

II, 1.63,0.693,0.210,0.OOOE+OO 
*HOURGLASS 

11,0,0,0,0,0,0 
*EOS_JWL 

11,3.71,3.230E-02,4.15,0.950,0.300,4.300E-02,1.0 

* SECTION_SOLID_ALE 

11 , 11,0 

0 

*PART 

material type # 8 (High Explosive Burn) 

11 , 11 , 11 , 11 , 11,0 

$ 

$ DEFINITION OF MATERIAL 12 (AIR) 

$ 

*MAT_NULL 

12,1.280E-03,0.000E+00,0.OOOE+OO,0.OOOE+OO, 0.OOOE+OO 
*EOS_LINEAR_POLYNOMIAL 

12,0.000E-00,0.0,0. OOOE+OO,0.OOOE+OO,0.400,0.400,0.OOOE+OO 
O.OOOE+00,0.OOOE+OO 
*SECTION_SOLID_ALE 
12,11,0 
0 

*PART 

material type # 9 (Fluid) 
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12,12,12,12,11,0 

$ 

$ NODES 

$ 

*NODE 

1,0.000000000E+00,0.OOOOOOOOOE+OO,0.000000000E+00,5,5 
2,0.000000000E+00,0.OOOOOOOOOE+OO,3.53 999996,5,5 
3,0.OOOOOOOOOE+OO,0.OOOOOOOOOE+OO,0.OOOOOOOOOE+OO,3,3 


11095,0.OOOOOOOOOE+OO,6.07999992,0.OOOOOOOOOE+OO,0,0 
11096,0.OOOOOOOOOE+OO,6.07999992,3.53999996,0,0 
$ 

$ ELEMENT CARDS FOR SOLID ELEMENTS 

$ 

* ELEMENT_SOLID 
1,11,1,125,127,3,2,126,128,4 
2,11,3,127,129,5,4,128,130,6 
3,11,5,129,131,7,6,130,132,8 


5367,12,10785,10909,10911,10787,10786,10910,10912,10788 

5368,12,10909,11033,11035,10911,10910,11034,11036,10912 

$ 

$ ELEMENT CARDS FOR SHELL ELEMENTS 

$ 

*ELEMENT_SHELL_THICKNESS 

1,1,11037,11039,11040,11038 

1.27,1.27,1.27,1.27 

2,1,11039,11041,11042,11040 

1.27,1.27,1.27,1.27 


28,1,11091,11093,11094,11092 

1.27,1.27,1.27,1.27 

29,1,11093,11095,11096,11094 

I. 27,1.27,1.27,1.27 

$ 

$ DETONATION POINTS 

$ 

*INITIAL_DETONATION 

II, 0.000E+00,0.OOOE+OO,0.OOOE+OO,0.000E+00 
*END 
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APPENDIX L. LS-DYNA INPUT DECK: 2-D AIR EXPLOSION WITH QUARTER- 

INCH THICK STEEL 


This Appendix provides the input deck used for the air explosion two-dimensional 
air explosion interacting with a half-inch thick steel. Appendix G and references 13 and 25 
offer the meaning to the command lines. 


*KEYWORD 
♦TITLE 

APM Steel-TNT Coupled Eulerian/Lagrangian 

* CONTROL_TERMINATION 

4000.,0,O.OOOE+O0,0.000E+00,0.000E+00 

* CONTROL_TIME_STEP 
0,0.67 

* DATABASE_BINARY_D3 PLOT 
50.0,0 

* CONTROL_ENERGY 

2 , 1 , 2,2 

* DATABAS E_BINARY_D3THDT 

1 . 0,0 

* DATABASE_HISTORY_NODE 
2000 

*DATABASE_NODOUT 

.8 

*DATABASE_GLSTAT 

1 . 

*CONTROL_ALE 
2 , 1 , 2,-1 
0. 

*CONSTRAINED_LAGRANGE_IN_SOLID 
1,2,0,0,0,4,1,0 
0 , 0 , 0.1 

*SET_PART_LIST 

1 

1 

* S ET_PART_L1ST 
2 

10,12 

*ALE_MULTI-MATERIAL_GROUP_PART 
10 

* ALE_MULTI -MATERIAL_GROUP_PART 
11 

*ALE_MULTI-MATERIAL_GROUP_PART 
12 

*ALE_INTERFACE_RECONSTRUCTION 

1 

$ 

$ DEFINITION OF MATERIAL 1 (Steel) 
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$ 

* MAT_PL AS TI C_jt INEMATIC 

1,7.87,2.1,0.300,4.300E-03,0.000E+00,0.000E+00 
0.000E+00,0.000E+00,0.0 
*HOURGLASS 

1,0,0.000E+00,0,0.000E+00,0.000E+00 
*SECTION_SHELL 

I, 2,0.000E+00,3.00,0.000E+00,0.000E+00,0 
0.635,0.635,0.635,0.635,0.000E+00 
*PART 

material type # 3 (Kinematic/Isotropic Elastic-Plastic) 

1 , 1 , 1 , 0 , 0,0 

$ 

$ DEFINITION OF MATERIAL 10 (WATER) 

$ 

*MAT_NULL 

10,1.00,0.000E+00,0.000E+00,0.000E+00,0.OOOE+OO 

* EO S_GRUNEISEN 

10,0.148,1.75,0.000E+00,0.000E+00,0.280,0.000E+00,0.000E+00 

1.00 

*SECTION_SOLID_ALE 

10 , 11,0 

0 

*PART 

material type # 9 (Fluid) 

10 , 10 , 10 , 10 , 11,0 

$ 

$ DEFINITION OF MATERIAL 11 (TNT) 

$ 

*MAT_HIGH_EXPLOSIVE_BURN 

II, 1.63,0.693,0.210,0.000E+00 
*HOURGLASS 

11 , 0 , 0 , 0 , 0 , 0,0 

*EOS__JWL 

11,3.71,3.230E-02,4.15,0.950,0.300,4.300E-02,1.0 
*SECTION_SOLID_ALE 
11 , 11,0 
0 

*PART 

material type # 8 (High Explosive Burn) 

11 , 11 , 11 , 11 , 11,0 

$ 

$ DEFINITION OF MATERIAL 12 (AIR) 

$ 

*MAT_NULL 

12,1.280E-03,0.OOOE+OO,0.000E+00,0.000E+00,0.000E+00 

* EOS_LINEAR_POLYNOMIAL 

12,0.000E-00,0.0,0.000E+00,0.000E+00,0.400,0.400,0.000E+00 
0.000E+00,0.OOOE+OO 
*SECTION_SOLID_ALE 
12 , 11,0 
0 

*PART 

material type # 9 (Fluid) 

12 , 12 , 12 , 12 , 11,0 
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$ 

$ NODES 

$ 

*NODE 

1,0.000000000E+00,O.OOOOOOOOOE+OO,O.OOOOOOOOOE+00,5,5 
2,0.000000000E+00,0.OOOOOOOOOE+OO,3.53999996,5,5 
3,0.000000000E+00,O.OOOOOOOOOE+OO,O.OOOOOOOOOE+OO,3,3 


11094,0.329162806,6.07108307,3.53999996,0,0 
11095,0.OOOOOOOOOE+OO,6.07999992,0.OOOOOOOOOE+OO, 0,0 
11096,0.OOOOOOOOOE+OO,6.07999992,3.53999996,0,0 
$ 

$ ELEMENT CARDS FOR SOLID ELEMENTS 

$ 

*ELEMENT_SOLID 

1,11,1,125,127,3,2,126,128,4 

2,11,3,127,129,5,4,128,130,6 

3,11,5,129,131,7,6,130,132,8 


5366,12,10661,10785,10787,10663,10662,10786,10788,10664 
5367,12,10785,10909,10911,10787,10786,10910,10912,10788 
5368,12,10909,11033,11035,10911,10910,11034,11036,10912 
$ 

$ ELEMENT CARDS FOR SHELL ELEMENTS 

$ 

*ELEMENT_SHELL_THICKNESS 

1,1,11037,11039,11040,11038 

0.635,0.635,0.635,0.635 

2,1,11039,11041,11042,11040 

0.635,0.635,0.635,0.635 

3,1,11041,11043,11044,11042 

0.635,0.635,0.635,0.635 


29,1,11093,11095,11096,11094 

0.635,0.635,0.635,0.635 

$ 

$ DETONATION POINTS 

$ 

*INITIAL_DETONATION 

11,0.OOOE+OO,0.000E+00,0.OOOE+OO,0.000E+00 
*END 
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